{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# _*Option Pricing with qGANs*_ \n",
    "\n",
    "\n",
    "## Introduction\n",
    "In this notebook, we discuss how a Quantum Machine Learning Algorithm, namely a quantum Generative Adversarial Network (qGAN), can facilitate the pricing of a European call option. More specifically, a qGAN can be trained such that a quantum circuit models the spot price of an asset underlying a European call option. The resulting model can then be integrated into a Quantum Amplitude Estimation based algorithm to evaluate the expected payoff - see [European Call Option Pricing](03_european_call_option_pricing.ipynb). <br/>\n",
    "For further details on learning and loading random distributions by training a qGAN please refer to <a href=\"https://www.nature.com/articles/s41534-019-0223-2\">Quantum Generative Adversarial Networks for Learning and Loading Random Distributions. Zoufal, Lucchi, Woerner. 2019.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:15.112855Z",
     "start_time": "2020-07-13T20:54:04.404067Z"
    }
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "from qiskit import Aer, QuantumRegister, QuantumCircuit\n",
    "from qiskit.circuit import ParameterVector\n",
    "from qiskit.circuit.library import TwoLocal, NormalDistribution\n",
    "from qiskit.quantum_info import Statevector\n",
    "\n",
    "from qiskit.aqua import aqua_globals, QuantumInstance\n",
    "from qiskit.aqua.algorithms import IterativeAmplitudeEstimation\n",
    "from qiskit.finance.applications import EuropeanCallExpectedValue"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Uncertainty Model\n",
    "\n",
    "The Black-Scholes model assumes that the spot price at maturity $S_T$ for a European call option is log-normally distributed. Thus, we can train a qGAN on samples from a log-normal distribution and use the result as an uncertainty model underlying the option.\n",
    "<br/>\n",
    "In the following, we construct a quantum circuit that loads the uncertainty model. The circuit output reads \n",
    "\n",
    "$$ \\big| g_{\\theta}\\rangle = \\sum_{j=0}^{2^n-1}\\sqrt{p_{\\theta}^{j}} \\big| j \\rangle , $$\n",
    "\n",
    "where the probabilities $p_{\\theta}^{j}$, for $j\\in \\left\\{0, \\ldots, {2^n-1} \\right\\}$, represent a model of the target distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:15.146155Z",
     "start_time": "2020-07-13T20:54:15.115128Z"
    }
   },
   "outputs": [],
   "source": [
    "# Set upper and lower data values\n",
    "bounds = np.array([0.,7.])\n",
    "# Set number of qubits used in the uncertainty model\n",
    "num_qubits = 3\n",
    "\n",
    "# Load the trained circuit parameters\n",
    "g_params = [0.29399714, 0.38853322, 0.9557694, 0.07245791, 6.02626428, 0.13537225]\n",
    "\n",
    "# Set an initial state for the generator circuit\n",
    "init_dist = NormalDistribution(num_qubits, mu=1., sigma=1., bounds=bounds)\n",
    "\n",
    "# construct the variational form\n",
    "var_form = TwoLocal(num_qubits, 'ry', 'cz', entanglement='circular', reps=1)\n",
    "\n",
    "# keep a list of the parameters so we can associate them to the list of numerical values\n",
    "# (otherwise we need a dictionary)\n",
    "theta = var_form.ordered_parameters  \n",
    "\n",
    "# compose the generator circuit, this is the circuit loading the uncertainty model\n",
    "g_circuit = init_dist.compose(var_form)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate Expected Payoff\n",
    "Now, the trained uncertainty model can be used to evaluate the expectation value of the option's payoff function with Quantum Amplitude Estimation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:15.153755Z",
     "start_time": "2020-07-13T20:54:15.149598Z"
    }
   },
   "outputs": [],
   "source": [
    "# set the strike price (should be within the low and the high value of the uncertainty)\n",
    "strike_price = 2\n",
    "\n",
    "# set the approximation scaling for the payoff function\n",
    "c_approx = 0.25\n",
    "\n",
    "# construct circuit for payoff function\n",
    "european_call_objective = EuropeanCallExpectedValue(\n",
    "    num_qubits,\n",
    "    strike_price=strike_price,\n",
    "    rescaling_factor=c_approx,\n",
    "    bounds=bounds\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot the probability distribution\n",
    "Next, we plot the trained probability distribution and, for reasons of comparison, also the target probability distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:15.764662Z",
     "start_time": "2020-07-13T20:54:15.156034Z"
    },
    "tags": [
     "nbsphinx-thumbnail"
    ]
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAERCAYAAABPbxE/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABQDklEQVR4nO2dd3yV9fX43ydhJGwIS2SKIFu2LAXEPaBOlLY42lK1tv21ar/WVkWt36/b2tZRR6VaFdyj1joZIkumrAAiyFJ2SAIJEO75/fF5Em6e3Jt7b3JHQs779XpeyfNZz3lubp7zfM7nfM4RVcUwDMMwEk1aqgUwDMMwagamcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykkHKFIyI9RORTETkgIttE5G4RSY/Qp6eI/Ndrf1BENonIsyJyXIi240RkuYgUisgqERmfuLsxDMMwwlErlRcXkabAJ8AqYBzQGXgYpwj/WE7XxsAG4AVgG9AJuBMYICKDVLXIG38E8AbwBPAr4DzgFRHZq6ofJeSmDMMwjJBIKjd+isjvgd8BHVQ11yv7HTAZaF1cFuVYZwIfAQNUdbFX9iFQW1VPD2r3H6CRqo4ob7zmzZtrx44dY7uhcti/fz/169eP23iJwuSMP9VFVpMz/lQXWeMp56JFi3apaouQlaqasgOYBUz1lbUHFLgwxrH6e/2GeOd1gUPAdb52E4EjQOPyxhswYIDGk+nTp8d1vERhcsaf6iKryRl/qous8ZQTWKhhnqupXsPpBmQHF6jqJuCAV1cuIpImInVE5CTgPuBLYIFX3Rmo7R8fWI0z2XWtnOiGYRhGLKRa4TQFckKU7/XqIvEf4CBOqTQDLlDVQNDYhBh/r6/eMAzDSAKpXsM5DNysqo/5yrcCU1T1DxH6d8Epmi44J4P9wHBVLRSR4cBsoK+qLvP1WQucpaof+8abBEwCaNWq1YCpU6dW9hZLyM/Pp0GDBnEbL1GYnPGnushqcsaf6iJrPOUcPXr0IlUdGKoupV5quNlGkxDljQk98ymFqq7zfp0vIp/jPNcmAP/g6EzGP37xeZnxVfVp4GmAgQMH6qhRoyKJEDUzZswgnuMlCpMz/lQXWaOR8/Dhw2zZsoXCwsLkCBWCxo0bk5GRkbLrx0J1kbUicmZkZNC2bVtq164ddZ9UK5xsfGs1ItIOqE/ZtZdyUdVvRWQPcIJXtB447I0/M6hpNyCAm+UYhhEDW7ZsoWHDhnTs2BERSYkMeXl5NGzYMCXXjpXqImuscqoqu3fvZsuWLXTq1Cnqfqlew/kAOFtEgu90PFBAaSUREc9xIAs3y0FVDwLTgct8TccDc1V1X0WFNoyaSmFhIVlZWRGVjSrMz4XLVkL9WZA2w/28fCUsyHX1RvVFRMjKyop5ppvqGc5TuA2Zb4rI/bjZyWTgEQ3agyMiXwMzVfUn3vlDQBEwH2ca647bz7MeCF54uQeYISJ/Bt7Gbfw8DzgngfdkGMc0kZTN4QBMzIZ3d0FhwJkTAA4E4I2d8J/dcGFzeKEb1E71K69RYSoyw03pn1tV9wJjgHTgPeAu4FFc1IBganltilkInAo8B7yPU1pv4Pbg7A8afzZwKXAG8CEwFpigFmXAMBKC6lFlcyBI2RQTAPYH4J1drl2sM52cnByeeeaZeIlbLm+//TarVq2Kqm3xgvu2bdu49NJLw7bLycnhiSeeKHesYcOGAW5N7YILLohSWodf5jvuuINPPvkkpjESSapnOKjqKuD0CG06+s6nUnomU17ft3GzGyMEp9+wqdT5FUMOcXdQ2WdPtE+2SEY1ZkEevOcpm/IoCLh2X+bB4EbRj5+Tk8Ozzz7Lb3/726j7lGw6TIvt/frtt9/mggsuoEePHlH3adOmDa+//nrY+mKFc8MNN5SpO3LkCOnp6cyZMycmOYPxy3z33XdXeKxEYBNawzBiRmaEPoYsdjOYaNgfgFMWhx4nHLfeeisbNmygb9++3HLLLeTn5zNmzBj69+9P7969eeeddwDYuHEj3bt354YbbqB///5s3ryZe+65h27dunHmmWdy5ZVX8tBDDwGwfv16zjnnHAYMGMCpp55KdnY2c+bM4d133+WWW26hb9++rF+/vpQcGzZsYOjQoQwaNIjbb7+9pHzjxo306tULgJUrVzJq1Cj69u1Lnz59WLduHbfeeivr168vkX/GjBmMHj2aCRMm0Lt3b4BS7sm5ublcdNFF9OjRg+uuu45AIFCmzeuvv87VV18dUuarr766RAF++umn9OvXj969e3Pttddy8OBBADp27Mi9995b8hlmZ8fkrxUTpnAMw6g23HfffXTq1ImlS5fy4IMPkpGRwVtvvcXixYuZPn06N910U3EIK9asWcPEiRNZsmQJO3fu5I033mDJkiW8+eabLFy4sGTMSZMm8de//pVFixbx0EMPccMNNzBs2DDGjh3Lgw8+yNKlS+ncuXMpOX79619z/fXX8+WXX9K6deuQsj711FNcf/31LF26lIULF9K2bVvuu+8+OnfuXCI/wIIFC7j33ntDmu8WLFjAww8/zPLly1m/fj1vvvlm2M+mPJkLCwu5+uqrmTZtGsuXL6eoqIgnn3yypD4rK4vFixdz/fXXlyjiRGAKxzCMaouqctttt9GnTx/OOOMMtm7dyvbt2wHo0KEDQ4YMAWD27NmMGzeOzMxMGjZsyIUXXgi4DY9z5szhsssuo2/fvvz85z/nu+++i3jdL774giuvvBKAH//4xyHbDB06lIcffpj777+fb7/9lszMzJDtBg8eHNa1ePDgwZxwwgmkp6dz5ZVXMnv27IiyhWLNmjV06tSJrl1dRK+rrrqKWbNmldSPHTsWgAEDBrBx48YKXSMaUr6GYxiGUVFeeukldu7cyaJFi6hduzYdO3YscdUNjn4cLqJKIBCgSZMmLF26NOZrR/LSmjBhAj179mTmzJmcffbZPPvss5xwwgll2pUXpdl/jeLz4PJoXJMjRZSpW7cuAOnp6RQVFUUcr6LYDMcwjJjRUaGPy1pE/1BJAy5vEXqccDRs2JD8/PyS83379tGyZUtq167N9OnT+fbbb0P2GzFiBO+99x6FhYXk5+fz/vvvA9CoUSM6derEa6+95u5LlWXLlpVcKy8vL+R4w4cPpzj01UsvvRSyzTfffEOnTp341a9+xdixY/nqq6/KHTMUCxYsYMOGDQQCAaZNm8aIES6rSqtWrVi9ejWBQIC33nqr1OcTavxu3bqxceNGvv76awBefPFFRo4cGbUc8cIUjmEYceOmdpAZ5VMlI821j4WsrCxOOeUUevXqxS233MIPf/hDFi5cyMCBA3nppZfo1i10kPlBgwYxduxYTj75ZC6++GIGDhxI48aNAacwnnvuOU4++WR69uxZ4nhwxRVX8OCDD9KvX78yTgOPPfYYjz/+OIMGDWLfvtB7yKdNm8Ypp5xC3759yc7OZuLEiWRlZTF8+PAS+SMxdOhQbr31Vnr16kWnTp246KKLALeWdcEFF3D66adz3HFHEx2HkzkjI4Pnn3+eyy67jN69e5OWlsZ1110X8frxJqXBO6syAwcO1OCFxcpSVeNpjb5hE7nN6rD5pEbsPi6DQLqQdkTJ2lZA+7V5fHlPa1IUwaRcqurnGYrqIms0cq5evZru3buHrVeFCavdPpuCcrzVMtNgXHN4uTsxf78qGi6mOEDlgQMHOO2003j66afp379/zOPEwrEa2qaYUN8HEamywTuNFHI4AKtOyWLX8ZkE0gTS3H9+oJaws209drfJZMJq2xFuRI+I+75MzHb7bAp8mz/TcDObsV6kgWS+zEyaNIlVq1ZRWFjIVVddlXBlY5TFFE4NpXhH+K7jMwnUCqFN0oRAmpTsCK/Im6hRM6md5r4vX+bBQ5tdKJuCgJvVnJ8FN7eDQTFs9owXL7/8cvIvapTCFE4NpXhHeEhlE0RFd4QbNRsR9315tWeqJTGqEmYoqaE8vLl8G3swBQHX3jAMozKYwqmhvL+7bGDFcAS89oZhGJXBFE4NJdrZTUXbG4Zh+DGFU0OJdq9EMRn2TTGqAJVJT3DeeeeRk5MTFzmCg2eGY/LkyfzlL38BIqcJiJQK4amnnuKFF14AYNSoUcSyZcOfEiFSCoVEYk4DNZTzs1wyrGgnLo1qQW6R+2kYxfjTW1SWSOkwyktPUBzePxz/+c9/Ki1fRYmUJqC8VAhFRUWV2qTpT4kQKYVCIrH31hpKLDvCAb4/BEMXw/qCxMlkGJHwpycIFd7/Bz/4AQMGDKBnz548/fTTJX07duzIrl27SlIX/OxnP6Nnz56cddZZFBS4L3aoVAUQPh2Bn3vvvZeTTjqJM844gzVr1pSUB6cJuPXWW+nRowd9+vTh5ptvDplWYNSoUdx2222MHDmSxx57jMmTJ5eK4vyvf/2LYcOG0atXLxYsWABQpk2vXr3YuHFjmZQIwSkUCgsLueaaaxgyZAj9+vVj+vTpAEyZMoWLL76Yc845hy5duvC73/2u4n+0IOx9tYYyuKFL8zvtuwCaHp3mWXUABi+CN3rCqKYJFtAwQnDffffx1VdflQTbnDFjBgsWLGDFihUlEZf/8Y9/0KxZMwoKChg0aBCXXHIJWVlZpcZZt24dr7zyCs888wyXX345b7zxBj/60Y+YNGkSTz31FF26dGH+/PnccMMNfPbZZyXpCCZOnMjjjz8eUrZFixYxdepUlixZQlFREf379y95sBezZ88e3nrrLbKzsxERcnJyaNKkCWPHjuWCCy4oZerKyclh5syZgFMmwezfv585c+Ywa9Ysrr32WlasWFHuZ7ZixYqSzyw4GnTxvcybN4+tW7dy1llnsXbtWgCWLl3KkiVLqFu3LieddBK//OUvadcuxlhEPmyGU0Mp3hFeLy9MZNiAklYUoKnvlWRPEZz5FTy9LfEyGkY0+MP7/+Uvf+Hkk09myJAhbN68mXXr1pXp06lTJ/r27QscDclfXqqCaNIRfP7551x00UXUq1ePRo0alYT8D6ZRo0ZkZGTw05/+lDfffJN69eqFva/x48eHrSuW5bTTTiM3N7fCa1OzZ88uuZ9u3brRoUOHEoUzZswYGjduTEZGBj169AgbGDUWTOHUYGqnQfph3yqOOkXTYssB+k3fwdahcEXL0k2KFH6+Fn61DorMe81IMcHh/WfMmMEnn3zC3LlzWbZsGf369QsZvr84HD8cDckfnKqg+Fi9enVJu0jpCKJpU6tWLRYsWMAll1zC22+/zTnnnBPVfUW6johQq1atkoygUPm0BaE+o8piCqcGk1sEeVl1S5Xdt20WI9/cQq95u2m09xCZ6S5MyT0dy/b/61Y4bznsPZwceQ3Dn57Az759+2jatCn16tUjOzubefPmRT12eakKoklHcNppp/HWW29RUFBAXl4e7733Xpk2+fn57Nu3j/POO48///nPJWauWNMWTJs2DXAzlMaNG9O4cWM6duzI4sWLAVi8eDEbNmyIOPZpp51Wcj9r165l06ZNnHTSSVHLESumcGowM3JA046+KdXLPUzTIwfLtBOBP3aE13tCPd835uO9Lo/92gOJldUwoGx6Aj/nnHMORUVF9OnTh9tvv70k42e0hEtVEE06gv79+zN+/Hj69u3LJZdcwqmnnlqmTV5eHhdccAF9+vRh5MiRPProo0D5qRBC0bRpU4YNG8Z1113Hc889B8All1zCnj176Nu3L08++WRJds/yUiLccMMNHDlyhCFDhjB+/HimTJlSamYTbyw9QRhqQnqCX66Dv209en78ujzuqDuXqfOO5hTxu6kuyYNxK2CzTy81qQWv9YAzmiVS4qNUxc8zHNVF1nikJ0gG1SXkP1QfWZOVnsBmODWYj/aUPm/2fWSf534NYUF/GOIL5JlTBOd8BY9vDd3PMAwj5QpHRHqIyKcickBEtonI3SISfveW6zNIRJ4Xka+9fmtE5E4RyfC1myIiGuIInRawBrGpENYG6RcJKE12ljWnhaJ1XZh+Mvy4VenyI8CN6+CGtS7XjmEYRjAp3YcjIk2BT4BVwDigM/AwThH+sZyu47229wPrgD7APd7PS3xts4FrfGUbKyl6tefjvaXPG+0+SK2i6M2rGenwz27Qqz7c+g0E93xyG2QfgNd6Qlbt+MhrGEb1J9UbP68DMoGLVTUX+FhEGgGTReQBrywU96vqzqDzGSJSCPxdRDqoarDD+H5Vjd5VpYbwsd+ctj2yC6UfEfhde+hWD364GvKPHK2bngOnLIL3ekP38N6dRjVEVaNyETaObSqy/p9qk9q5wIc+xTIVp4RGhuvkUzbFLPF+tgxRZwQRUPjEN8Np+n3sCqeYsc1hTj/omFG6fH2h82D7r6U2OGbIyMhg9+7dFXrYGMcOqsru3bvJyMiI3DiIVM9wugGfBReo6iYROeDVlXVkD88wXCzKNb7yHiKSC9QFvgT+oKozKy5y9WdJPuwO2sNV61CARnsPVWrM3g2cM8ElK+HzIK/R3CNw/nJ4uDP8uq2lqa7utG3bli1btrBzZ6h3vuRQWFgY84MuVVQXWSsiZ0ZGBm3bto2pT6oVTlMgJ0T5Xq8uKkSkNfAH4EXfbGkJMB+3RtQCuAlnthuhqgsqKnR1x29Oa7qjEInDC2uLOvDJyXD9WvjH90fLA8Bv1sOK/fBEV6iT6nm1UWFq165dKoxMKpgxYwb9+vVLqQzRUl1kTZacKd2HIyKHgZtV9TFf+VZgiqr+IYox6uAcD9oCA1R1bzltM3HKZ5mq/iBE/SRgEkCrVq0GFO8sjgf5+flR5dBIBr/lZJYE6fMJe1YxMn8LAM3qF7Jn/9E3na7t68Q8vgKv05an6EyA0lOaPuRwFytpQuXCE1SlzzMS1UVWkzP+VBdZ4ynn6NGjw+7DSfUMZy/QJER5Y0LPfEohbuXyBaAnMLw8ZQOgqgUi8h/gwjD1TwNPg9v4Gc/NelVl89+BI7ByNqXcyr6Z25jv9rsv2xVDsktv/JxYfn6ScIwGLtgNV6xyZrVivqIJv8kYzru9nBmuolSVzzMaqousJmf8qS6yJkvOVBs3snFrNSWISDugvlcXiUdx7tTjVDWa9sXU2BXPWTlwKOjuM/IPk7n/SNj2leHcLJjXHzr7TMMbC2HYEnhvV0IuaxhGFSXVCucD4GwRCY6pMB4oAMpd2BeR3wO/BH6kqrOjuZhnUjsXWFQxcas//v03FXGHjoXu9WH+ABjdpHR5/hEXIueBTWAOT4ZRM0i1wnkKOAi8KSJneGsok4FHghf/vYgCzwWdTwD+F2dO2yoiQ4KOFl6bxiLyuYj8XETGiMh4YDpwvNe3RuJXOJVxh46WrNrwYR+4rk3pcgX+5xu4OhsKEzPJMgyjCpHSNRxV3SsiY4C/4Vygc3Bmssm+prWA4HA3Z3k/r/aOYK4BpuAU2U5cxIKWQCEwFxipqvGLylmN+O4gLN9/9DwN56GWDGqnwRNdXGSCX69zYXCKeWE7rCuAN3u6sDmGYRybpNppAFVdBZweoU1H3/nVlFU0/j6FwMWVk+7Ywr/Zc1BDqH04efYsEfjF8dA1Ey5f5QJ+FjM3FwYvhnd7Qd+qH1zXMIwKkGqTmpFE/Oa0M5OUSsDPmc1gfn+neILZfBCGL4E3U7en0DCMBGIKp4agGkLhRL21Nv50rec82PwyHAi4aAV/2mjOBIZxrGEKp4awYj98HxS9pkF62Zw2yaZpbfhPb/jV8WXrbt8IE1ZDgTkTGMYxgymcGoJ/djOqSdUIMVMrDR7rAn/vCrV8cdam7oCRS2FbdGl6DMOo4lSBR46RDKqSOS0Uk9rAx32gmc+N5cs8GLQIFoZLVGEYRrXBFE4N4GAAZuaULqtqCgdgVFNYMAC61ytdvu0QjFji1nUuWwn1Z8HpjKT+LLh8JSzItfUew6gOmMKpAXyxDwqCUj4fX8clTauKdM6Euf3hPJ8H3UF16zpv7HSOBYpwIODOT1/q1nssrbVhVG1M4dQA/Oa0s5pV7bw0jWvBu73hphCpNvwTmQCwPwDv7IKJ2TbTMYyqTMwbP0WkNzAYaA1kAHuAtcCcSNGajdTgz39TFc1pftIFHjoRMtLg3k2R2xcEXDDQL/NgcIq97wzDCE1UCkdETgCuB34ItMK9WObgwsc0AeoBARGZCTwLTFNVM3BUAXYdgsX5pcvGVAOFU8zaAhCiC+9dEICHN8O0nomWyjCMihDRpCYizwIrgb7A3UA/IENVW6hqW1VtgItVdiGwHHgAWC0iIxImtRE1n+aUflj3bQAtY8+pljLe3x19LomA194wjKpJNDOcQqCbqn4broGq7sKlGvhARH4LXIaLymykGL857axqNLuB0s4OiWhvGEbyiKhwVPXGWAb0TGnTKiyRETeqWjibipCZ5rzSYmlvGEbVpFL/niLSS0R+ISI3es4ERhViXQFsCtqln5EGIxqnTp6KcH5WbF/S87MSJophGJWkwgpHRK4HZgGjgPOABSJyQ5zkMuLARz5z2qmNISM9dNuqyk3tYpu11K3C7t6GUdOJxmkg3BbB/wGGquplqnoe8AvgD/EUzqgcZfbfVDNzGsDghnBhc0gris6u9uIO+PPmBAtlGEaFiObdca2I/DBEueAcg4qxLXdViMMBmJ5TuixV+W8qgwi80A2aby1wSifg+5oFtMxuz9+sh6nbkyikYRhREY3CmQD8VkTmisigoPIHgHki8qqI/Bt4ArgvEUIasbMgD/KCQvu3rA2966dOnspQOw16zN9Nvxk7aLHlAGlFAUSVtKIALbYcoPuC3aT74tpMzIZPbRuyYVQpovFSmyUiA4GfAu+KyMfA/6jq4yIyGxjpNf2jqi5NnKhGLPjXb85oCmnVeH1DgEZ7DtFrnttoc8WQbKbO61ZSX6cwwMpRLSnyJjuHFS5aATP7Qj9LWW0YVYKolmPV8QxwErAdWC4itwHZqvoX71iaQDmNGAkVP+1Yptn2QqZ0K12WdwTO/Qo2FKRGJsMwShOTl5qq5qrqLcAQ4BQgW0QuTYhkRoXZV+RC9gdzRjV0GIiVH7aChzqXLtt+GM7+CnYeCt3HMIzkEZWXmoj8SUTmi8gSEXkaKFTVccDPgDtFZKaInJxwaY2omL4XgjMz96gHx9dNmThJ5aZ28FtflOl1BXD+csgvSo1MhmE4opnhPIeLk/YwcDsuSvTHIiKq+gkuxtprXtnTsQogIj1E5FMROSAi20TkbhEpd7eIiAwSkedF5Guv3xoRuVNEMkK0HSciy0WkUERWicj4WGWsbnxUzaMLVJYHO8OVLUuXfZkHl6+ynDmGkUqiUTjnAjer6quq+m/gKtxaTmcAVT2iqn/zymKylotIU+ATnEv1OFxw0JuAuyJ0He9d/37cptPHgd8CL/nGHwG8AUz37uN94BUROSsWOasbNW39xk+awJRuZc2IH+yBn62xnDmGkSqiCd6ZDfxYRBbhAnn+HNgPbAlu5OXC+XWM178OyAQuVtVc3CypETBZRB7wykJxv6ruDDqfISKFwN9FpENQoNHbgVmq+ivvfLqI9ATuAD6KUdZqwcYC+DpI7dcWGNkkZeKkjDpp8EZPGLUUlgSlZ/jndmhTF/73hJSJZhg1lmhmOFcBXYBdQB7OPfoyVS2Mw/XPBT70KZapOCU0MnQX8CmbYpZ4P1sCiEhdYDTwqq/dVGCoiFSzqGLR4Z/dDGsE9atZOJt40agW/Kc3dPIZWv9vE/x1S+g+hmEkjogKR1XXqOpQoCHQXFVPVNX/xun63XAzqODrbQIOeHWxMAwX+WCNd94ZqO0fH1iNu++usQpbHSizflPDzGl+WteFD/tA89qly3/9Nby2IzUyGUZNJWq3aFXdn4AU0k1xmUP97PXqokJEWuPiuL0YNFsq7u8ff6+v/pjhiJbdXV8d46fFmy713EynftC3XYEfrYYZFo3AMJKGaIQVVBH5MfCyqh4pt2HpPicCx6nq5xHaHcY5JDzmK98KTFHViMFARaQOzvGgLTCgWCmKyHBgNtBXVZcFte8CrAXOUtWPfWNNAiYBtGrVasDUqVMj3mu05Ofn06BBg7iNF4psGnI9A0rOG3KYt/iC8ixqazeV3qDSrH4he/YftUF1bV810oNGkhMiy7qAZtxGL44EvWfVp4jHWEJn9sdPWB/J+NvHA5Mz/lQXWeMp5+jRoxep6sBQddE4DdwE3CMiLwKvBz+8gxGRLOAc4ApcyoKfRDH2XqBJiPLGhJ75+K8pwAtAT2C4bwZW/Lt//OLzMuOr6tPA0wADBw7UUaNGRRIhambMmEE8xwvFnG+BDUfPz25RmzE9y7/m3TdsKnXuDxnz2cT2cZSw4kSSEyLLOgo47nsXZ62Y/dTi9jqDmNsfOpRxqo8PyfjbxwOTM/5UF1mTJWc0azh9cakIRgNLRCTX2wT6voi8KSKficgGYAfwGLAel5Lav1gfimx8azUi0g6oT9m1l1A8inOnHqeq/vbrgcP+8b3zAG6Wc0zhj59W0/bfRMOPW8P9Pg+17w7B2ctg9+HUyGQYNYVoY6lNU9UROG+1W4ClQBFOMWwH/omb3Rynqv9PVbdGef0PgLNFJDi84njcfp6Z5XUUkd8DvwR+pKqzQ8h8ELf/5jJf1Xhgrqrui1LGakF+EczxOZGbwgnNLe3g18eXLltTABcshwNRG44Nw4iVaExqJajqetzMIV48BfwKeFNE7gdOACYDjwS7SovI18BMVf2Jdz4B+F9gCrBVRIYEjbk+yG36HtwenT8Db+M2iZ6HU47HFLP2uQjJxZyYCZ0yUydPVUYEHjnRzWxeDXKwn5cL41fBWz2hVqWSrxuGEYqU/lt5ay5jgHTgPVyEgUeBO31Na3ltiimOFHA1MNd3nB80/mzgUuAM4ENgLDBBVY+5TZ/+/Tc2uymfNIEXusPoJqXL/70brltr0QgMIxHENMNJBKq6Cjg9QpuOvvOrccommvHfxs1ujmls/SZ26qbBW71g5BJYFuSk9tz3LhrB3Z1SJ5thHIuY4eAYYOtBWHXg6Hk6Zd/cjdA0rgUf9IEOvmja93wLT0a7EmkYRlSYwjkG+MRnThvcCJrUDt3WKMtxdeHDkyHLN9//xTp4M1QQJcMwKkRMCkdELhARU1JVjI/NnFZpTqoH7/eBer5oBBNWwec5qZLKMI4tYl3DeQfYISIv4CIBrE6ATEYMBNQcBuLFKY3g1Z4wbvnRBHYHFcaugM/7Qq8oNmKf7tugCnDFkEOlNq5+9kTV2ExrGMkm1tlKZ9xO/MuBFSIyV0R+5qUUMFLA8v2wI2jDYsN09+A0Ksb5WfDMSaXLcorgnK9gczzioxtGDSYmhaOqG1X1TlXtBJwJfI1zY/5ORF4UkdGJENIIj9+cNroJ1DajZ6W45ji41+ehtvUQnP0V7LFoBIZRYSr8aFLVz1T1x7gw/4uAHwKfiMgGEfmNiKTc5bomYOa0xPD79vCLNqXLVh+AscuhwKIRGEaFqLDCEZGRIjIFl3+mFy7N81nAa7gNnC/EQ0AjPIVHXISBYGp6/pt4IQKPdYFLW5Qu/yIXrlgFRYHUyGUY1ZlYvdQ6iMgdIrIe+Axohwvnf5yq/lJVP1XV3+GyhI6Lv7hGMLP3QWHQg69dXehq4WziRrrAi91gpC837Lu7ncu0RSMwjNiIdYbzDfAz4GXgRFUdo6qveIEyg1kJLIiHgEZ4/Oa0s5q6N3MjfmSkw9u9oHf90uVPfwd3f5samQyjuhKrwrkQ6KCqt6vqhnCNVHWtqpoDQYIps35j5rSE0KS2i0bQ3heNYPJGeHpbSkQyjGpJrApnINA6VIWIHCcid1ReJCMadhyCJflHzwUY0yRV0hz7HF8X/tsHmvlcYa5fC+/sSo1MhlHdiFXh3IlL5RyKNpSN8mwkiE99s5t+DaB51cgGfczSvT78uzdkBv3XBHBOBF8cU9mVDCMxxKpwBBfxIxRtOZrW2UgwZdZvzJyWFIY2hmk9Sv/jFAbgwuWwan/YboZhEEVoGxG5Cud1Bk7ZPCkivtySZAC9gWMuz0xVRC2cTUq5sDn8vSv8LChJ+V4vGkHbzHQybKOOYYQkms2ZB4Dd3u8C7AN8+9s5hEsX/UT8RDPCkX0AtgT5BWamwfDG4dsb8eenbVzG0Ds2Hi3bfBD2nNqCftO3U/uw+Uwbhp+ICkdVX8Nt5kREngfuLs9DzUg8/tnNaY1dMjEjufyxA2w7BE8Feartb1KH5SNacPLMHaTb5lDDKEVM4WdU9ZpECWJEj63fVA1E4G9d4PtD8HaQp9q+Fhl8dVpLah8MsPu4DKantyOtjZK1rYD2a/NQtf1SRs3E3ourGYcCMCOndJmt36SOdIGXu8MIn0kzp2UGO9tmEqiVBiIEaqWxs209loxqyYTVcNhmP0YNJBqngQXA1aq6SkS+JLyXGgCqOjhewhllmZcL+UFr0q3rQK/64dsbiSczHd7tBSMWw6qCoAr/NCZNCKQJ7+yCidlOUdlMx6hJRGNSWwkUBP1uq6EpxG9OO8PC2VQJmtaG/zsBxq2IbC8rCMB7u+DLPJcO3DBqCtE4DVwT9PvVCZXGiIg//81ZZk6rMvxrB+51LIoXgIIAPLwZpvVMtFSGUXWwNZxqxN7D7q04mDNM4VQZ3t8NpEU33QwUtzeMGkQ0azgR122CiXUNR0R6AH8FhgI5wLPAXaoadveciNQB7gWG4OK7Zahqmf90L1/PVf5yoLuqZsciZ1Xgsxz3oCqmV304rm641kayKYjRESDW9oZR3Yl2DSch6zYi0hT4BFiFy5/TGXgYN/P6Yzld6wE/xaVAmAOcXk7bbMDvzr2xYhKnFr85zbzTqhaZaXAgBiWSafYFo4YRzRrO1Qm8/nVAJnCxquYCH4tII2CyiDzglYWSKUdEmqmqisiNlK9w9qvqvPiLnnxC5b8xqg7nZ8Fr2zVqs1of8y40ahipfsc6F/jQp1im4pTQyPI6qtasfIvrC+CbwqPndQROa5IycYwQ3NQO0gLRfy3n5sHvv4EjNeqbbNRkIiocEVngrbMgIl9652GPGK/fDWfyKkFVN+Hit3WLcaxw9BCRXBE5KCKzRaRcRVZV8ZvThjeGeumpkcUIzeCG0HxrAWlF0dvV7tvkIk3nHE6gYIZRRUj1PpymOEcBP3u9usqyBJiPWyNqAdyEM9uNUNVqlQLbokNXfUSg+4LdrB6cxa7jMwmkSWnzWkARBdKl1D/RB3tg8GJ4p5fLuWMYxyqSSsuUiBwGblbVx3zlW4EpqvqHKMa4EfhrKC+1EG0zccpnmar+IET9JGASQKtWrQZMnTo1qvuIhvz8fBo0aFChvkeAHzCcfGqXlD3FQk4iP3ynKFm76VCp82b1C9mzP6PkvGv7qpHVLZKcUDVkXbvpEApsrNOIjxt2ZHlmcw5LOrX1CH0KdnJm3re0bV2HO+jJbkq7GNajiD+wmmGkxl+6Mt/RZFJd5ITqI2s85Rw9evQiVR0Yqi5mheO5JF8NDAaOA77DzSL+qaqHyukaaqwdwOOqepevPB/nGv1gFGNErXC89o8DF6pq+/LaDRw4UBcuXBjNkFExY8YMRo0aVaG+83NhyOKj51m1YMfwqNemy+X0GzaVOr9iSDZT5x21Zn72RLkfU9KIJCdUDVn9ckLoz3TbQbh4Bcz37asS4O6OcFuH+Px9Y6Ey39FkUl3khOojazzlFJGwCicmpwER6Q6sAx4HeuFevnt5518Xr/XEQDa+tRoRaQfUx7e2E2eq1TLtR771mzFNk/8wMuJLm7owoy9c07p0uQK3b4TLV0J+UQoEM4wEEquX2tO4BGydVXWIqo5V1SHAiV75UzGO9wFwtog0DCobj1szmhnjWBHxTGrnAoviPXYisfWbY5OMdHjuJPjrieD3/3hjFwxbAt8UhOxqGNWSWBXOQOAOz5OsBO/8DmBQjOM9BRwE3hSRM7w1lMnAI8Gu0iLytYg8F9xRRM4VkUuBvt75pd7RwTtvLCKfi8jPRWSMiIwHpgPHA/8bo5wpI68I5vp2I51p+W+OGUTgxrbw8cnOVBrM8v0waBF84s+vaxjVlFgVzkYgI0xdBlDWgF0OqroXGIN7wXsPuAt4FLjT17QWZV8Cn8RlIv2Jd16cmXS0d34Q2ImLWPAf3OwsBxipqvFbnEkwM3OgKMgA2DUTOoT7CxjVltFNYeEAONnnpbanCM7+Ch7dDDVr55lxLBJTxk/gVuBhEdmgqvOLC0VkCHA3cEusAqjqKsqPFICqdoymzFdfCFwcqzxVjY/MnFZj6JgJX/SHn2TDtJ1HywPAb9fDknz4e1eXf8cwqiMVCd7ZCJjjeZjtAFp6x27gNuDt+ItZcymzfmPmtGOa+unwSg/ouwlu21D6H+/F7bD6ALzVE9raLNeohlQkeOfKBMli+NhSCNkHjp6nA6ObpEoaI1mIwK0doE8DmLAK9gXFTV+YBwMXwRu9XLQJw6hOpDp4p1EO/tnNkEbQKFYjqFFtOS8LFgyAcStKv3hsPwyjl8LfusCkNikTzzBiJtXBO41yKLN+Y+a0GkfXejC/P1yYVbr8sMLP18L1a+GQ5dUxqgmmcKooAYVPzGHAwM1q3+4Ff+xQtu6pbTBmGWyPKcaHYaSGmBWOiIwXkU9EZJOI7PAfiRCyJrIsH3YFRRBulO6iERs1kzSBezrB6z2hvu+/dvY+t66zKC90X8OoKsQa2mYC8E/ga6At8C7wb2+cXOBv8RawpuJfvzm9KdSy+WiN55IWMLc/dPJ5qW05CCOWwEvbUyOXYURDrI+wW4B7gF9450+o6rVAJ2AXLo+NEQf88dPMnGYU07sBfDkAxjQpXV4YgB+thpu/hhhS8hhG0ohV4XQBvlDVI7jAnY0AVDUPuB+4Mb7i1UwKjjgzSTCmcIxgsmrDf/vAb9qWrXt4C5y3HPZYUjejihGrwtkHJUk8tgLdg+oEyCrTw4iZz/fBwaCdTx3qwomZqZPHqJrUSoNHToR/doO6vujhH+91cdhWVD5lkmHEjVgVzkKgj/f7u8AdIvIzEbkKeBCXF8eoJP71m7Oauc2AhhGKia3h835wvC//3DeFLo/SWztD9zOMZBOrwvk/jgbovANYADwBPI9bw5kUP9FqLrZ+Y8TKoEYu+OewRqXL9wfg4pVw5wbnam8YqSQmhaOq81R1mvd7jqqOAxoATVT1FFX9JhFC1iS2H4Kv9h89F5yHmmFEonVd+Kwv/Oy4snV3f+syjOZaUjcjhVRkH04dEZkkIs+KyPu4bJ/jvdTTRiXxb/Yc0NAtEBtGNNRNg6dPgie7QC2fGfad3c7Ets58SY0UkeoU04aPj33mtLNsdmNUgOuOh89Ohha+l5XVB2DwYvjv7tTIZdRsUp1i2ghC1fLfGPHj1CawaAD0b1C6PKcIzl8O938L8/bBZSuh/iw4nZHUnwWXr4QFuZbwzYg/qU4xbQSx6gB8FxQTq14aDLUQ9EYlaJcBs/vBD1uWLg8At25w0Qne3AkHAqAIBwLwxk44fSlMWA2HbQOpEUdSmmLaKI3fnDayibPJG0ZlyEyHF7vDQ53L/sMfwSmfYAI477Z3dsHEbJvpGPEj1sfZrcCfROSU4MKgFNP/Ey/BaiJl9t+YOc2IEyJwUzv4oA80iDJFdUEA3tsFX1pQUCNOWIrpKsLBAMzIKV1m+W+MeHNWMzi1EXywN3JbcErn4c0wrWdi5TJqBpZiuoowd5+zoxfTpg70qJc6eYxjl5n7IrcpJgC8bx5tRpywFNNVBL857YymFs7GSAwFMToCHDDHASNOVGhJWkTaiMglXhy1i0WkwpnVRaSHiHwqIgdEZJuI3C0i5VqZvc2nD4rI5yJSICJhlzVFZJyILBeRQhFZJSLjKyprIgkVP80wEkFmBf7rn9hqqayNyhPrxs90EXkC+BZ4Dfg78DrwrYg8LiKxjtcU+ARnshuHczy4CbgrQtd6wE9x+XfmlDP+COANYDpwLvA+8IqInBWLnIlm92FY6FuYPcMcBowEcX5WbP/4CvxiHZy0AKZ8Z7l2jIoT67vOXcC1OOeAjkCm9/M2r3xyjONd541xsap+rKpPedf4rYg0CtdJVXOAZqp6NvBWOePfDsxS1V+p6nRVvQX4L27PUJXhs72lF8n61IdWFijISBA3tavYLGdjIVyzBnovhNd2WDBQI3Zi/dpNBP6oqg+q6iZVPej9fBD3cL86xvHOBT5U1dygsqk4JTSyvI6q5e8OEJG6wGjgVV/VVGCoiFSZLZV+c5pFFzASyeCGcGHzyEqnlkAo23b2Abh8FQxY5BwKbJ+OES2xKpyWwFdh6r7y6mOhG5AdXOBFLTjg1VWGzkBt//jAatx9d63k+HFBtWw6Alu/MRKJCLzQDcY1h/ppZR8CabgoF5e2gG9OgV+0gdohHFiW5sMFy120ghlRulkbNZtYFc5a4IowdVcAa2IcrymQE6J8r1dXGYr7+8ff66tPKV8XwLcHj57XFTi1ysy9jGOV2mnwcneXzuCSFk7xCEp9T9HM6Auv9ID2mfC3rrB2MFzTOvQDY04ujF4GZy5zMdgMIxwSwTJVurHI5TiT1Gc4Z4HtuFnNZTjz1RWq+loM4x0GblbVx3zlW4EpqvqHKMa4EfirqoqvfDgwG+irqsuCyrvgFOdZqvqxr88kvCRyrVq1GjB16tRobyUi+fn5NGjQoEz527ThsaDJVn/28jDLyrRLFGs3HSp13qx+IXv2H41e1LV91VhMiiQnVA1Z/XJC1f1M/YT7jgaziUym0Inp5RgzhrGLa9lAZ/aHbVMZopGzqlBdZI2nnKNHj16kqgND1UWz8bMEVX1VRHJwC/uP4UxWh4FFwDn+B3gU7AWahChvTOiZT6xjE2L84vMy46vq07iI2AwcOFBHjRpVSRGOMmPGDEKN99gKXK5Uj8tPaMqo9vG7biTuvqF0+LsrhmQzdd5Ra+ZnE9snTZbyiCQnVA1Z/XJC1f1M/YT7jvqZCCzNg9s3wr9DbAqdQ3Pm0pwrWsJdHaFLnDcwRytnVaC6yJosOaM2qYlIbW/WsEJVh+IW9lsDmao6rALKBtz6Sqmnhoi0A+pTdu0lVtbjlKF/LagbbgP12kqOX2mKAs5DLRiLn2ZUB/o2hPd6w5x+cHqTsvUKvLIDui+An2bDpsJkS2hURWJZwzmCM6V1B1DVgKruUNXKeOV/AJwtIg2DysYDBcDMSoyLqh7E7b+5zFc1HpirqjEE+EgMC/Ig98jR8xa14eSqP/s2jBKGNoZP+8InJ8MpDcvWHwGe+x66zIdfrYPvD5ZtY9QcolY4nmJZB7SK4/WfAg4Cb4rIGd4aymTgkWBXaRH5WkSeC+4oIueKyKVAX+/8Uu/oENTsHmCUiPxZREaJyAPAebgNpinH7w49pimkWTgboxoypinM7Q/v9nL7yPwcUvjrVug8H37/Dew5nHwZjdQTq5faH4A7RKR3PC6uqnuBMTh3//dwa0OPAnf6mtai7JaAJ3HRDn7inb/mHaODxp8NXAqcAXwIjAUmqOpH8ZC/svjz39j+G6M6I+L29ywZCFN7QNfMsm0OBOC+TdBpHtyzEfKKki6mkUJichoA/ghkAUs9T7LtlN4kj6oOjmVAVV0FnB6hTcdoysL0fZsqmDIhtwjm+VxITeEYxwJpAuNbwiXN4YXtcNdG2OQzpeUegTs2wl+2wu/bw/VtXKI449gmVoWzEliRCEFqGtNznH27mG71XDpgwzhWqJUG1x4HP2wFz2yDP30L232mtF2H4ab1LufO7R1c+zqW5faYJVa36KsTJEeNw8xpRk2hbhrc2BauOQ7+thXu3wR7faa0bYfg+nXwwGaY3NEpqXRvPfP0Mi7xh0q5n3/2RNV0MzfKEpXCEZFM3GJ7R+A74FNV3Z5AuY55LH6aUdOonw7/0x6uawOPbIZHtkD+kdJtNhTCVdlOKd3dCS5u7mz2uc3qsPmkRuw+LoPp6e1Ia6NkbSug/do8VC13VHUhmhTTJ+BSCHQMKs4VkcuryuJ7dePbQlhbcPS8lsCoJikTxzCSSuNacFcn+OXxcP9mN+sp9G2uWHUALl0J/erDhlEtyG1Wl0CalLhxBmoJO9vWY3ebTCasdrHhapsprsoTzZ/oAdxGyVNxeWh6AktwuXCMCuA3pw1tBA1jXU0zjGpO8zrwYGdYf4pzGqgVYpayZD/ktMwkUCut7J6BNCFQK413dsHEbItaXR2IRuEMxaUk+EJVC1V1NfBzoL2IHJdY8Y5NzJxmGEdpUxee8AKEXtUq9r0aBQF4bxd8mRe5rZFaovnbHgd84ytbDwgutI0RAwGFT03hGEYZOmXClO6wYpCLWB0LBQHn6WZUbaJ9mbDJapxYkg+7gzx0GqfDwBAhQQyjptK9PrzWEzJimOoECB1I1KhaRLty8KGIhNoT/Km/XFVjTcJWo/AnWxvT1O1XMAyjNAdjjNJ4IAC3roerWjulZVQ9olE4dyVcihqErd8YRnRkpjklEgv3b3bHoIYwsRVc2QqyaidGPiN2IiocVTWFEycOHIEvfDGqz7R00oYRkvOz4I2dzlwWK1/mueO36904V7WC87IsikGqsY8/iczKcVFzi+mUAZ1DBDg0DANuaudmOZXhsMLbu+CildBmDvxyHSzMNRfqVGG7PxJEqHAc96zJhZMalZRZsjXDCM/ghi769KvbAm4fThgy0+D8Zm4m8+J2F6cwlD7ZXeQ2mf5tK/SoBxNbw49awfF1E3YLhg+b4SSRva1KR+c0c5phhEfERRBovrWAtKKA21MQTEBJKwowrjm83AOuPs4lg9s4BO7tFDo9QjGrDsCt30C7uXDWMnhpO+w/Er69ER9M4SSJfWl12N+kTsl5GqFT8xqGcZTaadBj/m76zdhBiy0HSCsKIOoUTYstB+g3fQev9Cgd1qZ9BtzWAbIHw9x+LopB0zC2HMU58vxoNbSeA9dmw4y9ZXWbER/MpJYkVmdklTof2BCamveMYUREgEZ7DtFrnttoc8WQbKbO6xa5n8CQxu549ES3T+ef38MHe6AohELJPwLPf++Ojhnw41bu6FIvzjdUg7EZTpJYlVla4Zxl5jTDSBp10+CSFvBub9g6FP58IvRrEL79xkK451vougCGL4a/b4McS4tdaUzhJAGl7AzH9t8YRmpoWQd+3RYWD4SvBsLN7aB1nfDt5+TCdWudye3ylfD+biiqiK+2YQonGexvXJvc9KOuMOmHAwxpVE4HwzCSQu8GLmL15iHwQW+4smX4kDoHFV7bCRcsh+Pnwm+/hmX5pduowvxcuGwl1J8FpzOS+rOcolpg7timcJLBHp93WpOdB20DmmFUIWqlwTlZztvt+2HwTFc4tXH49jsOw6NboO9COPlLl1BucwFMWA2nL4U3d7ooCYpwIOA2sJ6+1NUfrsGzI3vsJQgF9jWrw4qhzVnfp0mpunp5h2r8m45hVFUa14KftoFZ/Vyunskd4YSM8O2/2g83rYf28+HVHU7R+HVKANgfoMbn7jGFkwAOB2DVKVksHdWSncdnlkkctbVzwxr/pmMY1YETMuHOjvD1KfB5X/jpcdAoPXz7SP/SNT13jymcOKPq3mB2HR8mSyFYlkLDqGaIwIgm8MxJzuT2Snc4p1nFHqA1OXdPyhWOiPQQkU9F5ICIbBORu0WknHeIkn6NReR5EdkrIvtE5CURyfK1mSIiGuKI7MRfQRbkuTeY8kJxgL3pGEZ1JTMdrmgFH/SBLUPhwRPcXqFoCQCv74SbvoZ3d8HeGuRundKNnyLSFPgEWAWMAzoDD+MU4R8jdJ8GnAT8FPc3vB94GzjV1y4buMZXtrESYpfLw5udMomG4jedaT0TJY1hGInkuLpwc3v4nT8ncgQCwCNb3CFAn/owqgmMbAKnNTl2UyqkOtLAdUAmcLGq5gIfi0gjYLKIPOCVlUFEhgJnAyNVdZZXthWYLyJnqOonQc33q+q8xN7GUd7fHX049YDX3jCM6k1FcvcUo8Cy/e54bKsr61UfRjY+qoBalbNPqDqRapPaucCHPsUyFaeERkbot71Y2QCo6gJgg1eXMqKd3VS0vWEYVY/zs4hrALYV++HxbXD5KrfhtPsCuG4NTN0O3x2M22WSTqoVTjecyasEVd0EHPDqou7nsTpEvx4ikisiB0VktoiUp8gqTaz5Oyqb78MwjNRzUztIi1LhZKbB/3aCG4+H3lGmws4+AH//Dq5cDW3mQtf58LM18K/vYXNhJQRPMqk2qTUFckKU7/XqKtLvhKDzJcB83BpRC+AmnNluhDcjijuxZClM89obhlG9GdzQpVEo8U4NQ2YajGsOt7Z3nm8Auw/D5zkwcx/MzIGl+aHz+QSzrsAdz37nzjtlOPPbyMZuLahjhMSOqs7B6aHN8J/dUMBIMme559HN7VyKbonFEyJKRFPolysih4GbVfUxX/lWYIqq/iFMv4+BfFW9yFf+EtBRVYeH6ZeJUz7LVPUHIeonAZMAWrVqNWDq1Kkx39MqGnITfSkkoqMddTnCoyylO6lzVVu76VCp82b1C9mz/+gut67tq4bxOJKcUDVk9csJVfcz9ZOfn0+DBuVEtEwR1eU7unrTYZ7P6sWyzBYcljRUjioe0QC1VRmRtovfk02tclRKHrVYTmOW0YRlNGYdDQnE5AcHrSikDzmcTA592UcbCkpGKEL4P7oxh+YcRNAgQ5cQoC7KMCLLGY7Ro0cvUtWBoepSrXB2AI+r6l2+8nzgLlV9MEy/V4EWqjraV/4+gKqeX841HwcuVNX25ck2cOBAXbhwYXQ3EoSqC18RTZbCcc3h5e6JeZOIlrKZSUuHfv/siXI/pqQRSU6oGrL65YSq+5n6mTFjBqNGjUq1GGWoTt9RBfKa1WFT14bsbpOJpgtyRMnaVkD7NXksvLd1zOPmFsEX3uxn5j5YmBc6vUJ5tKnjOSA0dmkapueU7+RQmeeTiIRVOKk2qWXjW3MRkXZAfUKv0QT387s/4431dhTXTZiWLc5S+NmX3vQ6TUpv/gwoaQFlXJs0XuiWWmVjGEZ8qWjunvJoVAvOzXIHQH4RzM2FGTlOCS3Ig8MRnmjbDsErO9wRDcH7BAfHMdBwqpesPwDOFpGGQWXjgQJgZoR+rUVkRHGBiAzErd98EK6TZ1I7F1hUGaEjUZEshYZhGNHQoJZLT3/vCTC7P+SMgE9Phjs6uDWcunF6iU1ERIRUz3CeAn4FvCki9+MUxmTgkWBXaRH5Gpipqj8BUNW5IvIh8IKI3MzRjZ+zi/fgiEhj4N/Av4CvgebAb4DjgcsTfWOJeNMxDMPwUy8dTm/qDoDCI27WMzPHHXNyK7b9IhH7BFOqcFR1r4iMAf4GvIfzPHsUp3SCqQVlVuGv8Nr+AzdT+zdOeRVzENiJi1jQEigE5uI2i8a+OGMYhlENyEh3m0VPawK3A4cCzjQ2Mwf+sCG2seK9TzDVMxxUdRVweoQ2HUOU5eBC1vjD1hTXFwIXV15CwzCM6kudNBje2B33fhtbRIR47xO0VQTDMIwawvlZ0T/0E7FPMOUzHMMwUkNZd+ND3B1UVlXcjY34cVM7t9FzfxSznIw01z6e2AzHMAyjhjC4IVzYHNKKytc4mWkwtrmLOBBPTOEYhmHUEIr3CTbfWuCUjj/+W8Bt3xjXnITsEzSFYxiGUYNI5T5BW8MxDMOoYaRqn6DNcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykYArHMAzDSAqmcAzDMIykkHKFIyI9RORTETkgIttE5G4RSY+iX2MReV5E9orIPhF5SUSyQrQbJyLLRaRQRFaJyPjE3IlhGIZRHilVOCLSFPgEUGAccDdwE3BXFN2nAaOAnwJXA4OAt33jjwDeAKYD5wLvA6+IyFnxkN8wDMOInlopvv51QCZwsarmAh+LSCNgsog84JWVQUSGAmcDI1V1lle2FZgvImeo6ide09uBWar6K+98uoj0BO4APkrcbRmGYRh+Um1SOxf40KdYpuKU0MgI/bYXKxsAVV0AbPDqEJG6wGjgVV/fqcBQEWlcefENwzCMaEm1wukGZAcXqOom4IBXF3U/j9VB/ToDtUO0W427764VkNcwDMOoIKKqqbu4yGHgFlX9s698C/CCqt4Wpt/HwH5V/YGv/F/ACao6TESGA7OBfqq6NKjNicA64GxV/cjXfxIwyTs9CVhT8bsrQ3NgVxzHSxQmZ/ypLrKanPGnusgaTzk7qGqLUBWpXsMB5zDgR8KUV6Sf/1zC9VfVp4GnI1y3QojIQlUdmIix44nJGX+qi6wmZ/ypLrImS85Um9T2Ak1ClDcGcirQr0lQv71BZf42RBjfMAzDiDOpVjjZ+NZqRKQdUJ/QazRh+3kEr+2sBw6HaNcNCABrKyCvYRiGUUFSrXA+AM4WkYZBZeOBAmBmhH6tvX02AIjIQOAErw5VPYjbf3OZr+94YK6q7qu8+DGREFNdAjA54091kdXkjD/VRdakyJlqp4GmwCpgBXA/TmE8AvxZVf8Y1O5rYKaq/iSo7L84T7ObcTOW+4EdqnpqUJsRwAzgb7hNoed57c/xOwwYhmEYiSWlMxxV3QuMAdKB93ARBh4F7vQ1reW1CeYK3CzoH8ALwCLgIt/4s4FLgTOAD4GxwARTNoZhGMknpTMcwzAMo+aQ6jUcwzAMo4ZgCscwDEQkQ0SeFpEuqZbFOHYxk1qcEZEGuDhw3YCmuA2mOTh37Zmqmp866WJDRDKAll64oSqLiNQGjquKcopIJ6ATsEFVN6RYlnrlVDcBNuNiEc4GUNUDSRArZkSkI6Cq+m2qZQmHiIwGeuAcmpap6pwUi4SInAJ8qaqBlAmhqnbE4cBFMLgbyMd9yfJx/8BbgLygsrvwFH1VP4BLgCMpluEXuD1VecB84Mch2pySSjlxMfv+Auzx/sYPeOWPA0e8v/0R4CUgPYVyHolwBILPU/x3nwQ095X9GtgZJON24IYUy/kw8Keg8+OAud5nWQgc9GT9L9A4xbIGvM/sb8DwVMhQFULbHCtMBn6LUyhTVXVzcKWItMV51t2Jm/VMTrJ81Q4RuQL4K/AKsAQYDkwRkXE4xVOQSvmC+B0uL9MjOKXzSxFpAfwAl6tpMTACeAj4OfBESqR0+9tyPTl2++rq4x5ED1D+putk8SSwFC++lxfn8FFc9PfXvTaXAn8VkRxVfTkVQgKXA/8TdP43nNIZiTdT9H5/Aff9+AmpZSlwJXC9l9JlKu55tTgpV0+lxj2WDmArMCmKdpOArSmW9bMoj+WkduawEG+2EFQ2BvgemAdkeWWpnuFk44LQFp+PwL1N/j9fu7uAhSmUsw3wMu4h/kuCZlu4cFIB4LRUyeeTNQAM9n3GU0K0exFnJkqVnIXAqUHn+4HLQrSbAOyuCp8pbkY+1vsu5OJmYGtwL8MnJVIGcxqIH01wpp9IrCd0HLhkchrQCveWW96RlyoBPU4C/hNcoKqfAkNwD8i5ItI5FYL56AAsCDpf5P1c4Gs3G7e5OSWo6jZVnYDbr3YtsFxEzk6VPDHSGTfT9TMVt1aSKjbhvqfFFOEe4n5ygTpJkSgCqnpYVd/1vgutcMpwFXArsEpEFovILYm4tpnU4sc84HciMk9V94dqICL1cdPvuUmVrCwrgDWqOr68RiJyKS6Vd6rYhwubXgpV3Sgiw3Apw+cA9yRbMB/7Kf0ScdA7/Ivu6VSB/zlV/VxEBgDXAy+LyFzc+mNVIyPI0WE37k3czxEiR5ZPJP/AZSj+XFXX4ExnfxCR+aqaAyAizYDbKD9cV0pQZ5aeBkzzsi1fjDP9/wl4MN7XS/mX/xjiRuATYJOIfIgzAeTg/hma4LzWzsY9iMakRsQS5gPnRNFOOZrOIRUswq2DvO6vUNW9IjLGq/sLqX3orAEGAu94sgVwWWv99AQ2Jk+s8HgyPi4iU4F7gc9J7WcYiulBvwvOHPSJr00fnHNOqngQ93f9SkQ+wwUF7ot7Diz32vTBKcwfpUTCKFGXeXkKbp00KxHXMLfoOCIiTXBvjedw1C0aXKqEbFxg0aeK33xShWeG6qmq70Zol4lzi06J+6mIXAb8BrhAVfeEaZOOW2A+U1U7JVO+IBkuwa0nlRsAUUQ+ABaoqj90U8oRkR5AF+DzcJ91kuW5KkTxd1o2aeLrQLYGxV5MBSJyDnANbj2xNU5B7sVlGP438LSqptRELSLTgetVNWVOIaZwDMMwjKRgTgOGYRhGUjCFYxiGYSQFUziGYRhGUjCFYxiGYSQFUzhGWETkahFZJCJ5IrJXRJaIyCMJutblInJ1FO0mi4gGHdtE5I1oN4CKyBQRWVhpgeNAtPfstS2+73Vh6r/26icnSoYYxy31Ocf7OiKSJiI3et/JAhHJFZGVIvIXEamQK784loXxkCu+p45h6h4Xkecqct2ahCkcIyQi8nvgWVym1IuBibh9JmMTdMnLcXHHomEfMNQ7bsbte/jU21gbiXtiuE6iieWewYVR6SQiA4MLRWQQLtpBYRJkiBb/5xzv60zFbU58E/edvAoXlWKYVtz19nLcVoaKxGV7EPihiJxYwWvXCGzjpxGOG4G/q+ptQWXvichdqRIoiCJVnef9Pk9ENuE2Lp4HvOZv7O3VSVfVQ6oaTfihqsp+XCDQK3Bx5oq5Ahf7bkAqhComWZ+ziJwLXAacp6ofBFW9VdHZjcevgBdV9XDQtWrhlOePcbHorhSR9cBdqloShcOLfjEbtw/vpkrIcExjMxwjHE1wQTJLEfz2WGw2EZEfiEi2iBSKyGxvE2EpPJPKchE5KCKbReRe758ZEZmCS4UwMshUNjkGWYtjl3UMIddK3Jv/KcF1PtlOE5HpIpIvIvtEZIaI9AuqHyEiM0XkgIjsFpFnRKRheQKJyFARedcz+e0XkaUi8sPgz66C9zwVuLz4wer9vNwrj5sM3mfwum+8UV6bXsGfZaTPOdx1ROR8EQmIyxkUfJ1OXnm42fRI7+dn/oqKzm68mckwyka1+DUuGvhfcDOoa3HhbELtxH8DN8ux52oYbIZjhGMxLsz+JuDfquoPZ19MB1zY9dtx4e/vAj4UkS6qWgggImfh4jW9ANyCC/VxD+6f9jrv9/Y4JXeDN+6WGGTt6P383lf2AC5G2HYgZPIzERkFfIwLo3IVbhYxHDgeWCIiw4FPgbdx4fCzgPtwppdLy5GpA/AF8BTuQTwceF5EAqr6ChW/5zdxkRVG4GZ1pwItgLcoG/sqUTIE05HIn3O463wHbMN97pOD2l+Ny3tTKnBrEMWxCh8UkYfjFAljjDfuMl/5SOAzVX3Ae5H6QlU3hhljDi4YZu8Q4xhg6QnsCH3glMI3uPhaAWAl7qHSKKjNFK9+WFBZB1zE3OuCyuYB033j/w4XeLGtd/46MCMKuSbjwuvX8o6uOGWRi8v6GSxX3xD9pxCUIgAXSHUhYZLi4R7qftlP98bvFeVnKZ6sf8c9vIrLo7rn4Pv2fn8HeNz7/Qngbe/3XcDkeMgAzABe95WNCr7vGD/ncNf5E05JSZCcG4GHyvksWgNfeddWXDDa24AGlfi+P02INAfe57XZu+YUoGM5Y9Tyvvs/q6gcx/phUz8jJKr6FdAdtyD7BO5BcDuwUFwa7WJ2aFD6XHVvm4twgRaL7fr9Kbu2Mg1n0h1aAfGygMPesQYX8n+8qn4X1Garqi4tbxDPyeAU4J/qPTF89fU8+V4VkVrFBy7NwGHKWTMRkabiPKa+DZJ1Ek5BVpapwKUiUhc3yypjTkuCDMVE/Jwj8A/cS8oo73y0d/58uA6q+j3QDxcM90nczOleYI6I1IESD8ul3nHQM/kuFed1WTvEsK3xkr35uBc389mA+1+42Zv1hpKrCBewt3W5d1yDMYVjhEVVD6rqe6p6o6r2wGW17ELprIU7QnTdgct6CC69QG2cuSWY4vNmFRBtHzAIF6G5Le6t8wNfG//1QtEUp0i/K6c+HadwDwcdB3H31K6csacA43FmrrM8ef8BZEQhVyTeBRrgHob1gfdSIEMx0XzOYVHVb3CzqWu8omtwAU5XRuh3RFU/UtUbcOa653GmrKFe/RRV7Yt72SnCpVTuq6oDNMgpIIgM3N/Vf51N3rgX4Wb8I4DZEn57wEHi+/keU9gajhE1qvqciDyAi4RdTMsQTVviTHDg3hoPh2jXyvtZkcjERaoaaS9NNIvHe3HmwuPC1OdwNB14qPWEbaE6iUgGcD5wo6o+FVQelxc8Vd0vIv/GRdJ+TUPkX4qDDIWUTRgW6uUgHtF/nwWeEeeKfzExenmpakBEPsIpK//DvguwV8OvQRazhzAzE09B/VdcyvPJuJQZj4rInz2FFEwTKvadrhHYDMcIiYiUUSQi0gKXaTP4rbaluGRoxW3a494qF4B7E8WZ2C7zDXc57mFfnIzuEEl+M/Qe1POBicVeXyHq5+HS7i4McYRUOEBd3Myo5I3Z82rze11V5p6fxM1sngpTX1kZtlD6xQLgzApJWv51wDlCHMKZBtMIYyIEEJFWYarG4hLezfeVn0x0C/hrgDLpLUJ9L4AvvZ/NfG1bAPVwOXGMENgMxwjHchF5B/gIZyLrgNtkeQD4Z1C7XcCLIlLspXa3135KUJs7cZ5rz+MeJr1xnkvPqGqxV1Q2ME5EfoB72G0r54EeT27FJfX6QESextnrh+IWvP+Nc274VEQCuIXvPJwJ53zgD6pa5uGiqvtE5EvgDhHJxSnWW3GmwEZBTSt8z6o6A2eKCldfWRneAn4iIo/iMquOxq2ZVJSw96qqhSLyEvAL4BUtP1/UqyKSB7yKcy5oCfwQGIdbrPf3PRnnYBCJL3CfVQtV3RlU/rKILAFm4cyXA3Azy624XDfBDMTN+OZghCbVXgt2VM0D98//Ec5sVIj7534Z6BbUZgrOw+ti3FvdQdw/bhnvLdxawnLcm+wW3PpDraD65riH3B48M1YYuSbjeWuVI/sUgjykItXhXF9n4ZRpDs7rrW9Q/SnAf3GecPtx+d8fARqXI8OJuH0i+3F573/nlz3ae47hvkt5qVVWBuD3OA+tPOBfuFmE30stqs850r0CZ3jlZ0S4x2u9v8UW77u0B6cQR4Vp/x5wRRTf9zq4rJw/9pVf5F3ve5zSzsUp+n4hxngMn0ejHaUPS8BmVBhvQ18vVR0Yqa1hlIe3Njge6KQu/XW8xt0EnK2q/tlIqLaPASeq6vlh6qfgFOXGEHXpwLfArar6r0oJfQxjJjXDMFKGiJwE9MCFhLkrzsqmKW5TbLRrKg8Ca0Skq4YwlUbgMpxJOez6k2FOA4ZhpJa/40y1/8GFj4kbqrpXVTPVOa5E034LzuU/nNfi2ziTaygE+Im6vThGGMykZhiGYSQFm+EYhmEYScEUjmEYhpEUTOEYhmEYScEUjmEYhpEUTOEYhmEYScEUjmEYhpEUTOEYhmEYScEUjmEYhpEU/j+pk/ydeC2ECwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Evaluate trained probability distribution\n",
    "values = [bounds[0] + (bounds[1] - bounds[0]) * x / (2 ** num_qubits - 1) for x in range(2**num_qubits)]\n",
    "uncertainty_model = g_circuit.assign_parameters(dict(zip(theta, g_params)))\n",
    "amplitudes = Statevector.from_instruction(uncertainty_model).data\n",
    "\n",
    "x = np.array(values)\n",
    "y = np.abs(amplitudes) ** 2\n",
    "\n",
    "# Sample from target probability distribution\n",
    "N = 100000\n",
    "log_normal = np.random.lognormal(mean=1, sigma=1, size=N)\n",
    "log_normal = np.round(log_normal)\n",
    "log_normal = log_normal[log_normal <= 7]\n",
    "log_normal_samples = []\n",
    "for i in range(8):\n",
    "    log_normal_samples += [np.sum(log_normal==i)]\n",
    "log_normal_samples = np.array(log_normal_samples / sum(log_normal_samples))\n",
    "\n",
    "# Plot distributions\n",
    "plt.bar(x, y, width=0.2, label='trained distribution', color='royalblue')\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.grid()\n",
    "plt.xlabel('Spot Price at Maturity $S_T$ (\\$)', size=15)\n",
    "plt.ylabel('Probability ($\\%$)', size=15)\n",
    "plt.plot(log_normal_samples,'-o', color ='deepskyblue', label='target distribution', linewidth=4, markersize=12)\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate Expected Payoff\n",
    "Now, the trained uncertainty model can be used to evaluate the expectation value of the option's payoff function analytically and with Quantum Amplitude Estimation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:35.397951Z",
     "start_time": "2020-07-13T20:54:15.766834Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Analytically calculated expected payoff w.r.t. the target distribution:  1.0591\n",
      "Analytically calculated expected payoff w.r.t. the trained distribution: 0.9805\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEdCAYAAAAFP7AiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAogklEQVR4nO3deZyVZf3/8dcnlBBJTdHRFmfUr2bp128F2bdEgXEjt8KNctSvS46pmUsuKZWYkSWu5UqluIyNuaWigAtLKIiCWq78TAVE0ZRFHMYFmM/vj+s+eTjMdmbOuZc57+fjcR4z5z73ue8Px/H+nPu6rs91mbsjIiKV51NJByAiIslQAhARqVBKACIiFUoJQESkQikBiIhUKCUAEZEKpQQgZWdmo8zM8x5vmtmdZrZNQvGsb2aNZrY4iueoaPtxZvaama0ys6ltvHdIwb8l91gV578himW76LPdqGD7UVFM/eKOSbJlnaQDkIrxHjAs+n1r4ALgETPbwd1XxBzLCcD+wJHAG8ArZrY5cA1wJXA7sLSDY9QBr+Y9T6KgZjvgPGAcsCxv+/3At4Dm+EOSLFECkLiscvfHo98fN7MFwHRgH8IFN07bA3Pd/c7cBjMbBPQCrnf3f3biGP909+fKFWB3uPs7wDtJxyHppyYgScqc6GdN1CRzpZnNNbPmqBnmKjPbILezmd1uZlMKD2Jm55vZ22a2bvS8v5ndGDXvNJvZVDMbmLf/POBY4Gt5zTejCMkI4B/5zULFis53R8G2XLPRjtHzmuj5oWZ2nZm9Z2YLo3/Lpwreu5OZ3Wdmy8ysycyeMLM9zWwIcF+022vR8eZF71mrCaijzyX32ZjZxWZ2WhTP0qipbKOufBaSfkoAkpSa6OdbQF/Ct++RwHeAXwC1rHln8CdgsJltldtgZkZoxrnF3VdGm/8G7A2cAYwg/I1PMbP/il4fDjwAvERoJvlWdOyTotfrom33dxB/LzNbJ+/Rlf+XLgKagIOBW4BfRr/n/n3bA48BWwA/imK/G/gi8FT0bwQ4MIp5eDvn+hvtfy45hwK7A/XA2cB+wG+68G+TLHB3PfQo6wMYBbxLaHJch9B2PQVYDmzRyv7rALsQ2tW3jLZ9CpgPnJ+3X220z47R82HR88F5+6xPaA65Lm/bOGB2wTmH5B+rnX9Lbr/Cx6+j16cCd7R3bELyc+Cmgv2eARrznv8FWAis10Ys+0XHqSnYflS0vV+Rn8s84BVgnbxtlwNvJf03pEd5HuoDkLhsAqzMe74AGOHuiwDM7AjgdGBbwsUpZztggbu3mNk44EgzG+Xh6nQU4UKea4vfGXjH3afl3uzuK8xsPDCoxP+e7xMuljlvduEYDxY8fwHYMu95LeHu5oMuHDtfMZ/LFHfPH9H0ArCZmfV294+7GYekjBKAxOU9YA/CN9G3gDejizhmNhy4iTAK51xgCaHZ426gT94xbiA0Dw01syeBg/ikGYToPW+3cu63gY1L+Y8BnvfudwIvK3j+MWv+ezcBFnXzHFDc59JaTAb0jn6XHkQJQOKyyt1nt/HaIcAsdz8xt8HMBhfu5O7zzOxhwjf/rQjNQn/J22URsFkrx68iJJU4fEi4WObravJZTLh4d1caPhdJIXUCSxqsB3xUsK2ujX3/TPjmfyLwN3dflvfaLEJzxW65DWbWF9gXeLRk0bZvIWGYab49u3isR4BDzaxPG6/nvpG39XpOGj4XSSElAEmDh4DdzGykme1hZpcSRqK05m+Eb9lfJzQJ/Ye7TyKMmrnNzP7PzPYjjPhZDxhTruAL3A1sa2aXRf+W0YTRN11xPrAh8HczGxEd70wzOyZ6fW7083gz+6aZ/XdrB0nJ5yIppAQgaXAdcAlwCnAXUA0c1tqO7v4RMAF4HXi4lV2GExLK5YRhpAbUuvu/Sh516/HdT+jHOJiQDKqBU7t4rLmETtp3CUNV746OOz96fT6hD+RAwgX+vtaPBCT8uUg6WdQPJ5IJZrYO4QJ4vbv/Iul4RLJMncCSCWbWG/gfwp3BJoS7BhHpBiUAyYrPAU8A/waOd/eFCccjknlqAhIRqVDqBBYRqVBKACIiFSpTfQD9+/f3mpqaLr13xYoVrL/++h3vmAJZihWyFW+WYoVsxZulWCFb8XY31jlz5rzr7puu9ULSs9EV8xgwYIB31ZQpU7r83rhlKVb3bMWbpVjdsxVvlmJ1z1a83Y2Vgtlvcw81AYmIVCglABGRChVrAshbqq7w8aM44xARkeQ6gWuB/EUuXk0oDhGRipVUAnjS3ZsSOreIiKA+ABGR9GpogJoaBtfWQk1NeF5CSSWAV8xslZnNNbPjE4pBRCS9Ghqgvh7mz8fcYf788LyESSDuBLCIsKbrEcD+hJWKrjWz02KOQ0Qk3UaOhObmNbc1N4ftJZL4ZHBmdhthsfBN3b2lldfrgXqAqqqqAY2NjV06T1NTE/369etOqLHJUqyQrXizFCtkK94sxQrpj3dwbW345l/AzZg2eXJRxxo6dOgcdx+49sESru4lLAjuwNYd7atK4HTKUrxZitU9W/FmKVb3lMf7/vvuffq4w9qP6uqiD0cGKoE1L7WIyDvvQG0tfPQR9O695mt9+8Lo0SU7VRoSwEGENU/nJx2IiEii5s2DQYPg2Wfhnnvg+uuhuho3g+pqGDsW6upKdrpY6wDM7E7Cqk7/BHoBI6LHT7yV9n8RkYrx7LMwbFjo6H34Ydhll7C9ro5pU6cyZMiQkp8y7kKwucAxwBcBA14AjnT3m2OOQ0QkPR59FPbfPzTxTJ8OO+4Yy2ljTQDufi5wbpznFBFJtXvvhREjQhPPpEnhZ0zS0AcgIlKZ/vxnGD4cdtop3AXEePEHJQARkfi5w4UXwg9/CHvuCY88Av37xx6GEoCISJxaWuC00+Dcc+Gww0ITUEIFaUoAIiJx+fhjOPxwuOIKOPVUuPnmtcf6xyhTi8KLiGRWUxMcdBA8+CD89rdw1llglmhISgAiIuX2zjuw777w1FOhuOvoo5OOCFACEBEpr/nzYa+9YMECuPvuMN4/JZQARETKJb+696GHwjQPKaJOYBGRcnj0Udhtt/D79Ompu/iDEoCISOnde28Y319VBTNmxDa1Q7GUAERESun66+HAAxOr7i2GEoCISCnkqnuPPRb22COx6t5iKAGIiHRXiqp7i6EEICLSHR9/DEcckZrq3mJoGKiISFflV/deeCGcfXbi1b3FUAIQEemKXHXvnDlhWudjjkk6oqIpAYiIFKuwuveAA5KOqEuUAEREivHcc7D33qmt7i2GOoFFRDrr0Udh113D7ymt7i2GEoCISGfcd18mqnuLoQQgItKR668Pa/f+93+nvrq3GEoAIiJtcQ+Ltxx7LOy+O0yenPrq3mIoAYiItKalBU4/Hc45B37wg9AElIHq3mIoAYiIFMpV915+OZxyCtxyS2aqe4uhYaAiIvmamuDgg2HSpExW9xZDCUBEJOfdd0N17+zZma3uLYYSgIgIhOrevfcOPzNc3VuMRPsAzOzzZtZkZm5mPat3RUSy47nn4NvfhrffDhO7VcDFH5LvBB4DNCUcg4hUsvzq3r///ZPfK0BiCcDMdgWGARcnFYOIVKCGBqipYXBtLWy2GQwdGn7OmBEKvSpIIn0AZtYL+APwK2BZEjGISAVqaID6emhuxiBM6fypT4XVvHpIdW8xkroD+BHQB7gqofOLSCUaOTLM4pmvpSVU+1Ygc/d4T2i2CfAycLi7P2BmRwE3AJ9x97X6A8ysHqgHqKqqGtDY2Nil8zY1NdEvI1V8WYoVshVvlmKFbMWbhVgH19ZirVzz3IxpkycnEFHndPezHTp06Bx3H7jWC+4e6wO4FpiQ9/wowIF+Hb13wIAB3lVTpkzp8nvjlqVY3bMVb5Zidc9WvKmP9aOP3Pv2dQ8z/Kz5qK5OOrp2dfezBWZ7K9fUWPsAzGwH4BhgNzPbKNrcN/q5oZmtdvcP4oxJRCpArrq3uRnWXRdWrvzktb59YfTo5GJLUNx9ANsC6wIzgaXRI9cPsJDQMSwiUjrvvhtm8nzoIfjTn+CGG6C6GjcLHb9jx0JdXdJRJiLuUUCPAkMLtg0Dzgb2AV6NOR4R6cnaqu6tq2Pa1KkMGTIk0fCSFmsCcPd3gan528ysJvp1urfSCSwi0iXPPQfDhsGKFaG6t4IKvDor6UpgEZHSe+yxcMF3r7jq3mIkngDcfZy7m779i0hJ3Hcf7LFHqO597LGKq+4tRuIJQESkZG64Yc21e2tqko4o1ZQARCT73OF3vwvz9+fW7t1006SjSj0lABHJtpYW+OlP4Wc/67Fr95aLFoQRkez6+OPwrb+hAX7yE7jssjC5m3SKEoCIZFP+2r2/+U24A+iha/eWixKAiGRP/tq9f/oTHHts0hFlkhKAiGRLfnXvXXfBd7+bdESZpQQgItnx/PPh4t/UpOreElBviYhkw2OPwaBBYdTP9Om6+JeAEoCIpN/48Z9U91bg2r3logQgIuk2bhx873uw446q7i0xJQARSadcde/RR0NtLUyZoureElMCEJH0KazuHT9e1b1loFFAIpIuqu6NjRKAiKTHihVw0EGq7o2JEoCIpEN+de8f/wg//GHSEfV4SgAikrwFC2CvvVTdGzMlABFJlqp7E6OeFRFJzowZ4YLf0qK1exOgBCAiychV9266aUgEO+2UdEQVRwlAROKXq+7dYQdV9yZICUBE4uMOF130SXWv1u5NlBKAiMSjpQXOOAPOPhu+//3QBPSZzyQdVUXTKCARKb+VK0N17y23wMknw+WXq7o3BZQARKS8VqwIa/dOnAijR8M556i6NyWUAESkfBYvDtW9Tz6p6t4UivUezMwONrMZZrbYzD40s7lm9nMz6x1nHCISgwULwgpezzwDd96pi38KxX0HsAkwBRgDLAN2BkYBmwM/jjkWESmXwure3XZLOiJpRawJwN2vK9g0xcw2AE4ys5Pd3eOMR0TKYMYM2G8/6NMnVPeqwCu10tANvxhQE5BIVjU0QE0Ng2trw5q9Q4ZA//5hEXdd/FMtkU5gM+sFfBr4OvAT4Bp9+xfJoIYGqK+H5mYM4J13wvDO00+HrbZKOjrpQId3AGY22cy2j34/0sw2KcF5V0SP6cA04MwSHFNE4jZyJDQ3r7mtpQV++9tk4pGiWEdfvM1sJbCruz9uZquBb7n7E906qdnXgb6ETuBfAre6+4lt7FsP1ANUVVUNaGxs7NI5m5qa6JeRNUWzFCtkK94sxQrpj3dwbS3WyjXEzZg2eXICEXVe2j/bfN2NdejQoXPcfeBaL7h7uw/gVeASYEegBRgBfKWtR0fHa+X4RwIObNPRvgMGDPCumjJlSpffG7csxeqerXizFKt7yuP9+GP39dd3DzP8rPmork46ug6l+rMt0N1YgdneyjW1M30AFwJXA6dGF+pb29jPotd7dSYj5Xkq+rkV8EqR7xWRJOSqe1esgHXXDVM95PTtGyp+JfU6TADu/kczuxfYFvg7cBLwQglj2CX6+VoJjyki5VJY3bveejByJL5gAbblluHiX1eXdJTSCR0mADM7Erjf3R81s/OBe9z9za6czMwmAg8DzwOrCRf/nwK3ubu+/Yuk3YIFocDrtddCde/3vhe219UxbepUhgwZkmR0UqTONAHdAHyLMF7/l8AEoEsJAHgSOAqoAVYR+hfOAa7t4vFEJC7PPw/DhsH776u6t4foTAJYCnwu+j3Xzt8l7v4L4Bddfb+IJCRX3fvpT6u6twfpTAJ4GLjZzOZGz8eZ2Yq2dnb3nUsSmYikw/33wyGHwBe+AJMmqcCrB+lMAjgGOAHYnlC5+xrwTjmDEpGUuPFGOPZY+OpX4YEHwlQP0mN0ZhRQM6EOADPbAxjp7v8od2AikrAxY+Css2CPPeCuu7R8Yw9U1FxA7q57P5GerqUlXPgvuQRGjICbboLemq+xJyp6NlAz29rMrjGzZ83sjejn1Wa2dTkCFJEYrVwJRx0VLv4nnwy33qqLfw9W1B2AmQ0gLOjyITAeeBuoAg4C6sxsqLs/1c4hRCStVqwInb0TJsCvfw3nnqu1e3u4YqeDvhh4GvhO1DcAgJn1BR6IXq8tXXgiEov86t6xY+G445KOSGJQbALYGTg0/+IPoaPYzC4GbitZZCISj7aqe6XHKzYBfEBY17c1GxOahkQkK154IVz8Vd1bkYrtBL4f+K2ZDcrfGD2/ELivVIGJSJnNnAmDBsGqVaG6Vxf/ilNsAjidMH/PNDN7y8z+YWaLCKt6vUqY2E1E0u7++2H33WGTTcI0D5raoSIVWwewGBhkZsOAbwBbAIuAWe7+YBniE5FSu+kmOOYYVfdK0cNAe7n7anefCEwsU0wiUi6q7pU8xTYBvWFmF5nZl8sSjYiUR0sLnHFGuPiPGAHjx+viL0UngOuAg4HnzGyWmdWb2QZliEtESiW/uvfHPw7VvZ/+dNJRSQoUlQDc/Tx33xrYE5gLXAosMrOGaKI4EUmTFSvCuP6bbw7Vvb//PXyq6BlgpIfq0l+Cu0929yOBzYGTgS8Bk8xsnpmNMrPPtX8EESm7xYtDW//EiaG6d+RITe0ga+juV4GBwG6EtQKWAtOBHwL/MrPDu3lsEemq11+HXXeFp5+GO+7Q1A7Sqq7MBlptZueZ2SvAI4ShoMcAn3P3I4BqQl/BmJJGKiKd88IL8O1vwxtvhBW8hg9POiJJqWKHgU4mfONfCIwDbnD3+fn7uPtqM7sVOKVUQYpIJ82cGSZ1y63d+z//k3REkmLFzgX0LrAP8JC7t7c4/DOAFo8RiVNu7d7Pfz7M66O1e6UDxVYCH9rJ/VYC8zvcUURKQ9W90gXF3gEAYGZfALYD+hS+5u4PdDcoESnCxRfDmWeGuX3uvlsFXtJpxfYBfAb4K7BXblP0M785qFcJ4hKRjrS0wNlnhwQwYgTceKMKvKQoxY4CuhDYEtiVcPEfDgwB/gy8BvxvKYMTkTasXAlHHx0u/qrulS4qNgHsA4wGZkXP33T3v7t7PXAPcGYpgxORVuSqe2+6CS64QNW90mXF9gFUAa9HQz1XEFYBy3kAuLNkkYnI2hYvhv32gyee0Nq90m3Ffm14Hegf/f4ysF/ea9+kgyUhzewQM7vXzN4wsyYzm2NmPygyBpHKpOpeKbFi7wAeAvYA7gYuA240swHAR4QCsUs6eP/phL6C0/ikpuBWM+vv7n8oMhaRypFbu3f58lDdO3hw0hFJD9CpBGBm6xEu1v8GXjWzKne/2cyaCNNDrwf8mDAFRHv2d/d3855PjiaOOx1QAhDJaWiAkSMZvGBBGNPf1BSGd6q6V0qowwRgZlsDDwM1eZuXm9mh7n434W6gUwou/jlPA9/t7DFEeryGBqivh+bmMM767bfDLJ4XXKCLv5RUZ/oALgJaCEM/+wI7EC7aHX3b76xvAy+U6Fgi2TdyJDQ3r7nNHa64Ipl4pMey9qf0ATN7A/ipuzfmbdsOeBH4grsv6vLJzXYn9Csc4+7j2tinHqgHqKqqGtDY2Njabh1qamqiX79+XYw0XlmKFbIVbxZiHVxbi7Xy/6WbMW3y5AQi6pwsfLb5shRvd2MdOnToHHcfuNYL7t7ug/Dtf+eCbb2i7V/r6P3tHLcGeBu4u7PvGTBggHfVlClTuvzeuGUpVvdsxZv6WFevdt9gA/fwnX/NR3V10tG1K/WfbYEsxdvdWIHZ3so1tbPDQNu/TSiSmW0MTAAWAFo4RgQ+qe5dvhzWKeie69sXRo9OJi7psTo7DHSSma1qZfsjhdvdvd1pCM2sLzAe6A3s6+4rOhmDSM+1YgUcemiYyfOCC8JUziNH4gsWYFtuGS7+dXVJRyk9TGcSwPmlOpmZrQPcDmwL7OLu/y7VsUUyK7+697rrwggggLo6pk2dypAhQxINT3quDhOAu5csAQBXE+oJTgE2NrP8yeOedvePSngukfR7/fVQ4PXqq3D77XDggUlHJBWkS+sBdENuGunWxrNtBcyLLxSRhL34Iuy1l6p7JTGxJgB3r4nzfCKpNXNmaPbp3VvVvZIYzSErErcHHgird228MTz2mC7+khglAJE43XwzHHAAfPnL8OijsPXWSUckFUwJQCQul1wCRx4Z2vqnTIGqqqQjkgqnBCBSbu5w1llwxhlwyCGhCWiDDZKOSkQJQKSsctW9Y8bASSfBX/6itXslNeIeBipSOZqbQ3Xv/ffDr34FP/95mNZZJCWUAETKYcmSMMxz1qw1q3tFUkQJQKTUVN0rGaEEIFJK+dW9EyeC5vGRFFMCECmVxx+HffeFddeFadPgq19NOiKRdmkUkEgpTJgAtbXw2c/CjBm6+EsmKAGIdFeuunf77cPUDqrulYxQAhDpjlx17267wdSpqu6VTFECEOkKVfdKD6BOYJFirVwJxx0HN94IJ54Iv/899OqVdFQiRVMCECmGqnulB1ECEOms/Orea6+F449POiKRblECEOmMhQtDde+//qXqXukxlABEOvLii+Hi/957Ye1eVfdKD6EEINIeVfdKD6ZhoCJtmTAhrN2r6l7poZQARFpzyy2huvdLX1J1r/RYSgAihS69FI44QtW90uMpAYjk5Kp7f/pTOPhgVfdKj6cEIAJrrt174onQ2Ki1e6XHUwIQaW6G4cPD1A7nnw9XXqmpHaQixJ4AzOy/zOw6M/uHma02s6lxxyDyH0uWwJ57hhE/114Lv/ylpnaQipFEHcAOwD7A40DvBM4vEqi6VypcEgngPne/B8DM7gD6JxCDVDpV94rE3wTk7i1xn1OEhgaoqWFwbS1ssQV84xvw8cehulcXf6lQmgpCer6GBqivh+ZmDOCtt0I7/wUXqLpXKppGAUnPN3JkGOmTzx2uuCKZeERSwtw9uZNHfQDuPqSdfeqBeoCqqqoBjY2NXTpXU1MT/fr169J745alWCH98Q6urcVa+Tt3M6ZNnpxARJ2X9s82X5ZihWzF291Yhw4dOsfdB671grsn9gDuAKZ2dv8BAwZ4V02ZMqXL741blmJ1T3m8LS3uG2zgHr7zr/mork46ug6l+rMtkKVY3bMVb3djBWZ7K9dUNQFJz7VqFRxzDCxfDusUdHf17QujRycTl0hKKAFIz5Sr7h03LlT3jhsH1dW4GVRXw9ixUFeXdJQiiYp9FJCZ9SUUggF8HtjAzA6Onj/g7s2tv1Okk5Ysgf33h5kz4Zpr4Ec/Ctvr6pg2dSpDNOxTBEhmGOhmwO0F23LPtwLmxRqN9CwLF8KwYfDyy6G696CDko5IJLViTwDuPg/QZCtSei+9BHvtBcuWwcSJMHRo0hGJpJoKwaRnmDUrrN27zjqhuvdrX0s6IpHUUyewZN/EiVBbCxttFJZv1MVfpFOUACTbGhpCh29u7d5ttkk6IpHMUAKQ7LrsMjj8cNh1V63dK9IFSgCSPe5w9tlw+ulau1ekG9QJLNmyahUcd1wo7DrhBPjDH7R8o0gX6Q5AsiO/unfUKLjqKl38RbpBdwCSDUuWwAEHwIwZa1b3ikiXKQFI+qm6V6QslAAk3V56Kazdu3SpqntFSkwJQNJL1b0iZaVOYEmnXHXvhhuqulekTJQAJH1y1b3bbRc6fVXdK1IWSgCSLvnVvdOmqbpXpIyUACQd3OFnPwvVvQcdpOpekRioE1iSt2oV1NfDDTeE8f1XXqkCL5EY6A5AkpWr7r3hhlDde/XVuviLxER3AJKcpUtDZ++MGeHCf8IJSUckUlGUACQZb7wRCrxefhn++tcwq6eIxEoJQOKn6l6RVFACkHg98QTss4+qe0VSQJ3AEp/ct31V94qkghKAxCO/uldr94qkghKAlN/ll4fq3kGDwtq9m2+edEQighKAlFOuuve000J174QJoflHRFJBncBSHqruFUk93QFI6TU3w4EHhov/eeepulckpWJPAGb2FTN7xMyazexNM/uVmenq0FMsXQp77QXjx4cL/6hRYJZ0VCLSilibgMzss8DDwAvAd4FtgEsIiejnccYiZaDqXpFMifsO4EfAesCB7v6Qu18LnA+cbmblmfu3oQFqahhcWws1NeF5WmUpVlgz3s9/HnbaCRYsCJ29uviLpF7cCeA7wCR3X563rZGQFAaX/GwNDaEjcv58zB3mzw/P03hhzVKssHa8b74JS5bA2WeHpRxFJPXiHgW0PTA5f4O7LzCz5ui1+0p6tpEjQ4dkvuZmOPVU6N27pKfqtlNPzU6s0Hq8AH/8Y/jcRST1zN3jO5nZSuBMd7+8YPtC4CZ3P7eV99QD9QBVVVUDGhsbO32+wbW14dupxMbNmDZ5csc7JqSpqYl+/folHUanZSneLMUK2Yq3u7EOHTp0jrsPLNyeRB1Aa1dka2M77j4WGAswcOBAHzJkSOfPtOWWoSml0BZbwEMPdf44cdhzT1i0aO3taYwV2ozXttySov4bxWzq1Kmpjq9QluLNUqyQrXjLFWvcCWApsFEr2zcElpX8bKNHh3bq/KaKvn1hzBjYYYeSn65bxozJTqzQdryjRycXk4gUJe5O4JcIbf3/YWZfBNaPXiutujoYOxaqq3EzqK4Oz+vqSn6qbstSrJC9eEVkLXEngAnA3mb2mbxtI4APgGllOWNdHcybF9ql581L9wUqS7FC9uIVkTXEnQCuBT4C7jKzPaIO3lHApQVDQ0VEpMxi7QNw96VmtjtwJWHI5zLgMkISEBGRGMU+CsjdXwBUKSQikjDNBioiUqGUAEREKlSslcDdZWbvAK1UdnVKf+DdEoZTTlmKFbIVb5ZihWzFm6VYIVvxdjfWanfftHBjphJAd5jZ7NZKodMoS7FCtuLNUqyQrXizFCtkK95yxaomIBGRCqUEICJSoSopAYxNOoAiZClWyFa8WYoVshVvlmKFbMVbllgrpg9ARETWVEl3ACIikkcJQESkQikBiIhUKCUAkZQysz5mNtbMtk06FumZelwnsJn1AwYTFp75LGGpyWWEBWemuXtTctF1jZn1ATZz9wVJx9IRM1sX2CKtsZrZVsBWwGvu/loK4unbzssbAa8D3wEeBXD35nb2T5SZ1QDu7l2t1o+FmQ0FvgK0AP9w9xkJhwSAmX0TeNLdW2I7qbv3iAdhXeFfAU2E/7BNhP95FgLv5207nyjxZeUBHASsTkEcJwGvRJ/nLOCIVvb5ZtKxAusCvweWRP/NL4q2XwWsjv4WVgMNQK+EY13dwaMl/3kK/gbqgf4F204B3smL823gxBTEegnw67znWwAzo8/0Q8LaJKuBicCGKYi3JfrsrgR2ieOcSSwKXy6jgNMJF/hGd389/0Uz+wLwfeA8wl3BqJjjyzQz+z7wB+AvwNPALsA4M/suIRF8kGR8Bc4CfghcSkgCJ5vZpsD3gKOAp4BBwMXA8cDViUQZfAAsj2JZXPDa+oSLwUWUY8nUrrkGeIZoXppoUafLgL8Cd0T7HAz8wcyWufutSQQZORQ4O+/5lYQkMJjojir6/SbC38qxsUbXumeAHwAnmNkbQCPhevZUWc6WdNYrYfZ8A6jvxH71wBtJxxvFMrmTj2dJ/lv1bKJv0nnbdgfeAh4HNom2peEO4CXgzLzngwjfrk4t2O98YHbCsX4OuJVwQT2ZvDsSYMMo7t2SjLEg3hZg54LPelwr+91MaM5IMtYPgV3znq8ADmllv8OAxWn5bAl3sAdEfxfLCXcpcwlfXr9UynP2pE7gjQjNEx15Jdo3DXYDqgjf/Np7vJ9UgHm+BDyQv8HdHwH+l3Chmmlm2yQRWCuqgSfyns+Jfj5RsN+jwNaxRNQGd3/T3Q8DhgPHAM+a2d5JxlSkbQh3hYUaCe3sSVpA+LvNWUW4oBZaDvSOJaJOcPeV7n5v9HdRRUhQLwA/A14ws6fM7MxSnKsnNQE9DpxlZo+7+4rWdjCz9Qm3hDNjjaxtzwFz3X1EezuZ2cHAbfGE1Kb3CFPSrsHd55nZt4H7gRnABXEH1ooVrJnkP4oehR2ovUjJ/wPuPt3MBgAnALea2UxCn1Ya9cnrvF5M+IZaaDWhqTVJ1wOjzGy6u88lNPWMNLNZ7r4MwMw2Bs4FpiUXZts8NK3eBtxmZhsABxKasn8NjOnu8VPxx18iPwYeBhaY2STCrekywh/hRoRRQXsTLgS7JxPiWmYBwzqxnxM6uZM0h9CGfkfhC/7JWs93EDpfk/4ffy4wELgHwMOoivVa2W8HYF58YbUvivMqM2sERgPTSf6zbM2UvN+N0GzxcME+OxEGYSRpDOG/8T/NbDLw/4CvEq4Rz0b77ERIYocnEmER3H05MI7Q97ZJKY7Zo4aBmtlGhG9Qw/hkGCjAUkJCmABcm8v+SYuaTHZw93s72G89wjDQxIbXmdkhwGnAfu6+pI19ehE6Cfd0963ijK8gjoMIfRLtTqBlZhOAJ9z9vHgiK46ZfQXYFpje1mceNzP7v1Y2L3L3Bwv2uwN4yd1/Hk9kbTOzYcDRhP6pzQlJaynwIjAeGOvuiTezmtkU4AR3j63Dv0clABER6bye1AksIiJFUAIQEalQSgCSSWZ2lJnNMbP3zWypmT1tZpeW6VzbmdmoqI+po31HmZnnPd40szs7M0TWzMaZ2eySBC3SCUoAkjlmdg7wJ2ASYVjckYQRPweU6ZTbEYpwNurk/u8B34oeZxBGnjwSDUNuzwWESmWRWPSkYaBSOX4MXOfu5+Ztu8/Mzk8qoAKr3P3x6PfHzWwBYUjnPsDthTub2Xru/oG7d6aQUaRkdAcgWbQRYQqKNXjekDYzq4maYA4zs5ujpqJ/m9laQz7NrNbMZpnZh2b2tpldHc0qi5kNAe6Ldn0tOua8IuPNVSLXRMecZ2aXmNkvzGwhUXVqa01AZlZtZn8xs3fNrNnM/mlmh+W93sfMLjKz183sIzP7h5ntU2R8UqF0ByBZ9BRhgrcFwHh3L5xELd8YwljvgwlTb5xnZu+6+1Xwn7H2E4GHCLOufhH4LWGKiGHRuc4gTNZ2ILCIUExYjJroZ37SOgx4HjiRNv4/NLPNCFXrzVEMrwM7RjHm3EEoxDqPMM3JocC9ZjbQ3Z8pMk6pMEoAkkUnAX8jVEW6mb0I3AlcHFVL5nve3Y+Pfp8UXVTPNbNrosrbXwLzgQPcfTWAmS0hlN5/y91nmtnc6P1Pu/u8zgRoZrn/t7YmzDb6PmtXy+7n7h+2c5jTCPMsDXD3RdG2R/LOsTuwLzDE3XNTGTxoZtsBI4FDOhOrVC41AUnmuPs/gS8TOn2vJlR2/gKYnWu6yXN3wfO7CDNwfiF6vjNwd+7iH7mTMHHYoC6GuAmwMnrMJSSBEXkXcYBHOrj4A9QCEwvel28Pwl3FY2a2Tu5BSBIDuxi7VBDdAUgmuftHhLb5+wDM7FjCyKBjgSvydv13wVtzz7cgzBa5BWERjvxjrzazxcDGXQzvPcLF2QkX6Dd97ZL7t9d619o2AZ5s5/X+hKkNVrbyWmsTtImsQQlAegR3/7OZXUSYAyrfZm08X5T3c419ojmNNiEsJtMVq9y9o/H8nZmDZTEhQbVlCWEdjO91Mi6RNagJSDInascv3LYpob288Jv18ILnuY7chdHzWcDw6KKfv886fLJq1MfRzz7dCLsrHgH2NrOqdl7fHGhy99mFj/jClKzSHYBk0bNmdg/wIKFJp5owSqYZuLFg3x3M7DpCu/5uhCaiU/yThbd/TVji8m9mdg2hb+B3wCR3z60bkesEPj6aqrnZ3Z+l/C4jFLlNN7PRhFFAXwbWd/eLCCOXJgEPmdnvCKOKNiAUnvVx93NiiFEyTAlAsuhXwHcJaw9sTGhnn0HoaH2tYN+zgP0ICeBDQrXtlbkX3f15M/sO8BtCB/FywgpXZ+XtM9/MzgB+Qli2cSGfDO0sG3d/x8x2IawJfDnwaeBl4MLodTezAwkLmpwKbEloFnqGsH6zSLs0HbT0SGZWA7wG7O/u4xMORySV1AcgIlKhlABERCqUmoBERCqU7gBERCqUEoCISIVSAhARqVBKACIiFUoJQESkQikBiIhUqP8PnYy6SulUXxoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Evaluate payoff for different distributions\n",
    "payoff = np.array([0,0,0,1,2,3,4,5])\n",
    "ep = np.dot(log_normal_samples, payoff)\n",
    "print(\"Analytically calculated expected payoff w.r.t. the target distribution:  %.4f\" % ep)\n",
    "ep_trained = np.dot(y, payoff)\n",
    "print(\"Analytically calculated expected payoff w.r.t. the trained distribution: %.4f\" % ep_trained)\n",
    "\n",
    "# Plot exact payoff function (evaluated on the grid of the trained uncertainty model) \n",
    "x = np.array(values)\n",
    "y_strike = np.maximum(0, x - strike_price)\n",
    "plt.plot(x, y_strike, 'ro-')\n",
    "plt.grid()\n",
    "plt.title('Payoff Function', size=15)\n",
    "plt.xlabel('Spot Price', size=15)\n",
    "plt.ylabel('Payoff', size=15)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# construct A operator for QAE\n",
    "european_call = european_call_objective.compose(uncertainty_model, front=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set target precision and confidence level\n",
    "epsilon = 0.01\n",
    "alpha = 0.05\n",
    "\n",
    "# construct amplitude estimation \n",
    "ae = IterativeAmplitudeEstimation(epsilon=epsilon, alpha=alpha,\n",
    "                                  state_preparation=european_call,\n",
    "                                  objective_qubits=[num_qubits],\n",
    "                                  post_processing=european_call_objective.post_processing)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "result = ae.run(quantum_instance=Aer.get_backend('qasm_simulator'), shots=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact value:        \t0.9805\n",
      "Estimated value:    \t1.0070\n",
      "Confidence interval:\t[0.9806, 1.0334]\n"
     ]
    }
   ],
   "source": [
    "conf_int = np.array(result['confidence_interval'])\n",
    "print('Exact value:        \\t%.4f' % ep_trained)\n",
    "print('Estimated value:    \\t%.4f' % (result['estimation']))\n",
    "print('Confidence interval:\\t[%.4f, %.4f]' % tuple(conf_int))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-13T20:54:35.972442Z",
     "start_time": "2020-07-13T20:54:35.400448Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<h3>Version Information</h3><table><tr><th>Qiskit Software</th><th>Version</th></tr><tr><td>Qiskit</td><td>0.23.3</td></tr><tr><td>Terra</td><td>0.16.2</td></tr><tr><td>Aer</td><td>0.7.3</td></tr><tr><td>Ignis</td><td>0.5.1</td></tr><tr><td>Aqua</td><td>0.8.1</td></tr><tr><td>IBM Q Provider</td><td>0.11.1</td></tr><tr><th>System information</th></tr><tr><td>Python</td><td>3.8.5 (default, Sep  4 2020, 07:30:14) \n",
       "[GCC 7.3.0]</td></tr><tr><td>OS</td><td>Linux</td></tr><tr><td>CPUs</td><td>2</td></tr><tr><td>Memory (Gb)</td><td>3.7363624572753906</td></tr><tr><td colspan='2'>Thu Jan 28 21:29:04 2021 IST</td></tr></table>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<div style='width: 100%; background-color:#d5d9e0;padding-left: 10px; padding-bottom: 10px; padding-right: 10px; padding-top: 5px'><h3>This code is a part of Qiskit</h3><p>&copy; Copyright IBM 2017, 2021.</p><p>This code is licensed under the Apache License, Version 2.0. You may<br>obtain a copy of this license in the LICENSE.txt file in the root directory<br> of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.<p>Any modifications or derivative works of this code must retain this<br>copyright notice, and modified files need to carry a notice indicating<br>that they have been altered from the originals.</p></div>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import qiskit.tools.jupyter\n",
    "%qiskit_version_table\n",
    "%qiskit_copyright"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
